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ABSTRACT 

Composite  reliability  is  strongly  dependent  on  the  fiber  strength 
distribution.  Current  methods  of  gathering  statistics  through  single  fiber  failure 
methods  are  inefficient  and  costly.  This  thesis  develops  a  testing  method  from 
which  the  fiber  statistics  in  the  form  of  Weibull  parameters  can  be  accurately 
extracted  from  bundle  failure  tests. 

The  values  obtained  from  the  bundle  experiment  as  compared  to  known 
single  fiber  test  bench  mark  parameters  were  practically  indistinguishable.  The 
confidence  of  the  results  stems  from  a  thorough  analysis  of  the  bundle 
mechanisms  and  minimized  contaminates  which  can  disturb  the  strength 
distribution. 
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I.    INTRODUCTION 

A  composite,  by  definition,  is  any  two  or  more  materials  combined  on  a 
macroscopic  scale  to  form  a  useful  material  [Ref  1].  Today,  in  one  form  or 
another,  composites  are  being  used  on  every  level  of  our  society.  This  thesis 
will  focus  on  a  class  of  composites  known  as  "fibrous  composites,"  today  these 
are  typically  used  in  high  tech  applications,  often  high  performance  military 
aircraft.  Fiber  composites  consist  of  very  thin  fibers  with  high  aspect  ratios  and 
high  strengths  imbedded  in  a  matrix,  commonly  a  polymeric  matrix  such  as 
epoxy.  Glass/epoxy  is  an  example  of  a  frequently  utilized  composite  mixture 
that  has  been  in  use  for  many  years. 

One  of  the  many  advantages  this  class  of  material  provides  is  high 
strength-to-weight,  though  due  to  a  lack  of  thorough  understanding  in  failure 
mechanics,  composites  have  been  limited  to  secondary  load  carrying  structural 
designs  such  as  skins  and  control  surfaces.  In  recent  years  (last  20),  composites 
are  expanding  into  a  much  wider  role. 

A.    REQUIREMENTS  IN  ADVANCED  STRUCTURES 

Since  the  introduction  of  glass/epoxy,  advances  in  production  techniques  and 
developments  in  the  textile  industry  have  resulted  in  fiber  materials  with  much 
improved  strengths  and  stiffness.  Two  of  the  early  break-throughs  were 
graphite  and  boron  fibers  which  have  typical  strengths  and  stiffness  much 
higher    than    convention  ductile    isotropic  metals.      These    fibers  allow  for  new 


design  limits  but  also  give  rise  to  the  requirement  of  new  design  technique 
utilizing  a  material  of  an  unconventional  material  redundancy  with  anisotropy. 

Employing  the  fiber/matrix  as  unidirectional  lamina  to  form  multidirectional 
laminate  gave  the  designer  the  tools  to  make  structural  designs  subjected  to 
multiple  loading  conditions  a  manageable  task.  This  idea  shifted  the  attention 
of  the  designer  from  the  microscopic  level  to  the  macroscopic  level. 

Utilizing  the  high  strength  and  stiffness,  fibers  application  for  structural 
enhancement  were  incorporated  in  many  of  the  military  designs  in  the  early 
70's  such  as  the  F-lll  aft  fuselage  structure.  Testing  of  these  designs  were 
accomplished  by  component.  This  method  is  time  consuming  and  expensive  but 
must  be  done.  To  gain  any  type  of  reliability  assurance  many  tests  are 
required.  It  soon  becomes  apparent  that  testing  costs  could  out  pace  production 
costs. 

B.    FUTURE  DESIGN  REQUIREMENTS 

As  weaponry  becomes  more  and  more  sophisticated  the  demands  on 
composite  materials  for  a  wide  variety  of  uses  are  ever  increasing.  New  designs 
often  are  only  limited  by  the  structural  limitations  of  the  materials  in  use. 
Such  restrictions  might  include  the  thickness  of  a  wing  or  stabilizer,  the  sweep 
of  the  wing,  or  the  "g"  loading.  Lately  new  design  specifications  have  included 
the  use  of  composites  for  the  radiation  energy  absorbing  properties.  Though  the 
Navy's  interest  in  composites  spans  a  multitude  of  uses,  one  special  interest  is 
in  the  area  of  structural  enhancement,  such  as  rocket  motor  cases,  pressure 
vessels  for  submarine  flasks  and  jet  aircraft  pilot  ejection  seat. 


In  the  last  few  decades  composites  have  advanced  dramatically,  and  are 
beginning  to  provide  improved  alternatives  to  more  conventional  methods  of 
design  normally  utilizing  isotropic  alloys.  Fiber  reinforced  plastics  have 
attracted  a  large  amount  of  attention.  Improved  techniques  in  the  textile 
industry  have  lead  to  more  uniform  fibers  from  bases  of  precursors  for  high 
performance  graphite's  with  strengths  and  stiffness  exceeding  high  yield  steels. 

As  utilization  goes  up  the  increased  use  of  new  fibers  in  composites  and  the 
cost  of  system  failure  testing  demand  better  methods  for  reliability  evaluation. 


II.  BACKGROUND 

A.    COMPOSITE  LOAD  SHARING 

A  brittle  material  generally  has  a  higher  strength  to  density  ratio  but  a 
lack  of  ductility  limits  the  structural  uses.  Ductility  resists  crack  propagation 
from  a  region  of  damage.  Metals  have  a  high  ductility  and  do  provide  practical 
strength  and  stiffness  properties  for  many  structural  applications.  Fibers  do  not 
use  ductility  to  resist  failure  [Ref.  2:pp.  1-19],  because  ductility  is  not  required 
the  fiber  can  be  synthesized  from  brittle  materials. 

The  maximum  theoretical  elastic  strain  for  carbon  hydrogen  bonds  is 
approximately  emax  =0.1,  which  is  one  order  of  magnitude  greater  than 
present  day  strong  elastic  reinforcing  fibers. 

Presently  fibers  can  not  obtain  strengths  any  where  near  the  theoretical 
ideal  strength  of  the  material  but  are  an  order  of  magnitude  greater  than 
conventional  bulk  materials.  Fiber  strength  can  also  be  attributed  to  the  small 
fiber  diameter,  which  limits  the  flaw  sizes. 

Composites  are  geometrically  advantageous  to  homogeneous  material  such  as 
an  all  metal  or  all  ceramic  material,  for  its  ability  to  utilize  the  high  strengths 
and  stiffness  of  brittle  materials  while  maintaining  strength  redundancy  through 
the  matrix. 

The  strength  interaction  of  the  fiber  and  matrix  is  the  key  to  reliability 
effectiveness  of  a  composite.  Consider  several  fibers  in  parallel  without  matrix. 
The  load  is  equally  shared  by  each  fiber.  If  one  fiber  should  fail  the  failed  fiber 
is   completely   ineffective   and   no  longer  supports   any  load.      The   load  of  the 


remaining  fibers  increases  equally.     With  a  matrix  present  (Figure  1)  a  crucial 

phenomenon  occurs.     At  the  point  of  fracture  the  more  ductile  matrix  absorbs 

the  load  through  shear  stresses  along  adjacent  fibers  and  back  to  the  broken 

fiber  itself.      Near  the  fracture  site  high  shear  stress   develope  in  the   matrix 

along    the    broken    fiber    but    quickly 

dissipate      as      the      normal      stress 

increases.      On    the   two    surrounding 

fibers   a  noticeable   stress  increase  is 

observed     adjacent    to    the     fracture. 

The    short   distance    from    the   broken 

end    of   the    fiber    to    the    point    that 

fiber  is  carrying  a  full  stress  is  called 

the  ineffective  length  [Ref.  3] 

The    ability    of   the    composite    to 

adjust  to  multiple  fiber  failures  is  it's 

strength      redundancy,      and      is    Figure  1.        Local    load    sharing    in    a 

Fiber/Matrix  composite, 
dependent    on    the    characteristics    of 

both  the  matrix  and  fiber.  There  are  several  fracture  mechanisms  all  of  which 
are  a  function  of  the  bond  strength  (shear  strength  of  interface),  the  matrix 
shear  strength  and  the  distribution  of  the  fiber  strength.  The  most  critical  is 
the  fiber  strength  distribution.  The  distribution  of  the  fiber  determines  the 
composite  reliability. 
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B.    STATISTICS  AND  RELIABILITY 

The  complexities  of  composites  make  the  task  of  analysis  virtually 
impossible.  New  fibers  are  being  developed  every  year  and  attempting  to 
generate  thorough  material  testing  studies  as  has  been  done  with  conventional 
materials  is  not  realistic.  But  composite  reliability  can  be  estimated  through 
probability  studies.  Probability  of  the  composite  must  be  based  or  inferred  from 
accurate  statistics  of  the  fiber  distribution. 

To  determine  the  reliability  of  the  composite,  the  strength  distribution  of 
the  fiber  must  be  known.  The  distribution  can  be  quantified  by  strength 
statistics  (a  Empirical  Probability  Distribution  Function  (EPDF),  or  Empirical 
Cumulative  Distribution  Function  (ECDF)).  For  computational  efficiency  the 
empirical  data  can  then  be  modeled  with  an  analytical  function  and  thereby 
characterize  the  distribution  by  parameters  of  the  function.  To  define  the 
parameters  the  distribution  function  must  be  identified  by  model  type  (i.e., 
Normal,  Weibull,  etc.).  Using  the  model  Probability  plot;  the  closeness  of  fit  of 
the  data  can  be  judge  and  therefore  the  appropriateness  of  the  model.  Selecting 
the  model  allows  for  the  calculation  of  the  model  parameters. 

The  relationship  between  the  fiber  statistics  and  the  composite  probability  is 
represented  in  a  linearized  Probability  plot  (Figure  2). 

Reliability  is  mostly  concerned  with  the  lower  tail  of  the  distribution,  this 
also  being  the  most  difficult  to  define.  Because  of  the  redundancy  provided  by 
the  matrix  induced  load  sharing  the  composite  is  considerably  more  reliable 
than  it's  constituent  fibers  in  the  lower  tail.  As  the  fiber  statistics  fluctuate  in 
the  extreme  lower  tail  the  composite  reliability  fluctuates  with  it.    This  can  be 
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Figure  2.     Comparison  of  Fiber  and  Composite  in  a 
Probability  Plot. 


shown  in  a  representation  of  a  probability  plot  in  Figure  3.    There  is  a  critical 
point  where  the  fiber  plot  stabilizes,  the  linear  region  in  this  plot  determines 
the    accurate    statistics    of    the    lower    bound    reliability    (worse    case)    of    the 
composite. 

Fibers  have  a  characteristic  which  helps  to  identify  the  analytical 
distribution  model  that  will  closely  define  the  empirical  statistic.  Traditionally, 
the  weakest  link  model  [Ref.  4]  depicted  the  fiber  as  being  made  of  many  small 
segments  linked  together  much  like  a  chain.  The  segments  each  have  an 
intrinsic  strength,  none  being  exactly  equal.  When  the  fiber  is  stressed  the 
weakest    will    fail  and  the    total    fiber    fails.      The    Weibull  distribution  model 


defines    this    weakest   link   behavior   and   is   commonly   used   in    modeling   fiber 
strength. 

1.     Weibull  Density  Distribution  Function 

The  Two  parameter  Weibull  cumulative  distribution  function  (CDF)  is 
defined  by  the  following  equation. 

F(x)  =  1  -  exp  {-(x/fi)°),    x,  B,  a  >0  (2.1) 

where:  x  =    random  variable  (stress,  load) 

13  =  scale 
a  =  shape 

The    parameter    alpha    defines    the    shape    (spread    and    skew)    of   the 

distribution.  Alpha  less  than  3.5  the 
plot  is  skewed  positively,  greater  than 
3.5  it  is  skewed  negatively  and 
neutral  at  3.5.  The  alpha  parameter 
determines  the  strength  scatter  of  the 
fiber,  an  increasing  alpha  defines  less 
scatter  and  higher  strength 
uniformity.  The     beta     parameter 

defines  the  central  tendency  of  the 
distribution  and  is  relatable  to  the 
mean.  Alpha  for  a  fiber  type  is 
constant    for    the    weak    link    physics, 
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figure  3.    Extreme  lower  tail. 


whereas  beta  changes   with  the   fiber  length.      An  increase  in  length  causes  a 
decrease  in  beta  (a  weaker  fiber). 


The  composite  reliability  can  be  inferred  from  the  distribution  of  a  fiber 
strength.  To  gain  high  reliability  with  high  confidence  levels  a  proportionate 
amount  of  statistics  must  be  obtained  (i.e.  for  a  reliability  level  of  103  requires 
104  statistics). 

Presently  fiber  statistics  are  gathered  by  single  filament  testing  (Figure 
4).  First  a  single  fiber  specimen  is  prepared,  since  a  fiber  has  such  a  small 
diameter,  (AS-4  diameter  -  10  micron)  it  is  very  fragile.  The  fiber  is  placed  in 
a  cardboard  frame.  Once  the  specimen  is  loaded  in  the  test  apparatus  the 
cardboard  frame  is  severed  by  a 
heated      wire      followed      by      load 


fiber 


application   and   interpretation   of  the 

load    deformation    data.       This    is    a 

time      consuming      effort      which 

produces  only     one     failure     strength 

statistic.        To    gain    the    number    of 

figure  4. 

statistics      necessary      would      take     specimen. 


Single  fiber  testing 


r 

\ 
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case 
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Single       filament       test 


thousands    of   tests,    as    an    example    man-safe    reliability    is    106    therefore    to 
achieve  the  necessary  statistics  would  require  107  tests. 

A  new  approach  being  explored  is  to  test  many  fibers  in  parallel  (a 
bundle)  and  recover  the  fiber  characteristics  from  the  failure  test.  This  would 
provide  a  large  statistical  base  required  and  could  be  accumulated  with  a  much 
higher  efficiency.  With  a  validated  method  the  analysis  of  fiber  reliability  could 
keep  pace  with  the  development  of  new  fibers  with  reduced  effort. 


2.     Size  Effect. 

The  extreme  lower  tail  is  difficult  to  determine  accurately  for  the  fiber. 
For  a  Weibull  failure  physics  and  the  resulting  strength  distribution  the  lower 
tail  can  be  estimated  by  assuming  weakest  link  [Ref.  41.  The  parameter  beta 
is  a  function  of  the  fiber  length  and  is  defined  by  equation  2.2 

J32  =  ^(lAJ™-  (2.2) 

where  13,  =  beta  for  1, 
fl2  =  beta  for  12 
a    =  a,  =  a2 

A  detailed  development  can  be  seen  in  Appendix  A.  Using  this 
relationship  parameters  for  beta  can  be  retrieved  analytically  where  practical 
testing  is  impossible  and  the  extreme  lower  tail  is  estimated. 

This  thesis  will  attempt  to  measure  the  lower  tail  via  bundle  testing,  as 
opposed  to  the  aforementioned  estimation,  (fiber  statistics  via  single  filament 
testing).  This  investigation  explores  the  methodology  needed  in  testing  a  bundle 
of  fibers  to  failure.  By  recording  the  load  and  displacement  it  is  expected  that 
the  underlying  distribution  of  the  fibers  in  that  bundle  may  be  extracted.  If, 
for  example,  a  bundle  has  3k  of  fibers,  the  distribution  is  based  on  3k  of 
statistics.  As  the  statistical  base  approaches  the  complement  of  the  reliability, 
the  accuracy  of  the  estimation  increases. 
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III.  TECHNICAL  ISSUES  AND  RESOLUTIONS 

A  bundle  of  fibers  tested  to  failure  presents  several  experimental  obstacles 
which  effects  the  analysis.  There  are  three  experimental  consequences  to 
testing  a  bundle  in  parallel  that  are  not  experienced  in  testing  a  single  fiber  in 
series.  These  phenomenon  which  may  disturb  the  failure  distribution,  are  1)  a 
protective/adhesive  enhancement  coating  called  sizing  which  is  placed  on  the 
filament  spools  during  production,  2)  friction  which  occurs  during  testing  after  a 
portion  of  the  bundle  has  failed  and  3)  slack  which  results  from  specimen 
fabrication.    Each  of  these  is  discussed  in  the  following  paragraphs. 

A.  SIZING 

Fibers  are  packaged  on  spools  of  bundles,  in  order  to  prevent  the  fibers 
from  entangling  a  small  amount  of  a  chemical  liquid  called  sizing  is  coated  on 
the  bundle.  The  effect  of  this  sizing  on  bundle  testing  is  unknown  and  must  be 
address  at  sometime.  The  testing  in  this  experiment  was  done  leaving  sizing  on 
the  bundle.  A  method  to  remove  the  sizing  chemical  without  damage  to  the 
fibers  needs  to  be  studied  but  is  not  addressed  here. 

B.  FRICTION 

During  the  bundle  test  prior  to  fiber  failure  there  is  negligible  friction,  this 
is  due  to  the  equal  strain  between  fibers.  As  the  fibers  fail  the  strain  of  the 
failed  filaments  return  to  zero  leaving  the  loose  ends  entangled  in  the 
remaining  fibers,  this  places  added  strain  on  those  fibers.  The  resulting  effect 
shows  up  as  reduced  failure  strengths  and  clustering  of  the  these  failures.     In 
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the  relatively  short  gauge  lengths  friction  is  expected  to  be  small  and  very 
pronounced  in  the  long  gauge  lengths.  The  effects  of  lubricating  the  bundle 
with  oil  during  testing  will  be  observed  to  see  if  the  effect  is  reduced. 

C.    SLACK 

During  specimen  preparation  and  placement  in  the  test  apparatus  the 
alignment  of  the  fibers  becomes  misaligned  with  unequal  length.  The  effective 
result  is  an  initial  nonlinear  load/displacement  which  also  effect  the  failure 
distribution.  Another  concern  is  the  deformation  in  the  actual  load  train  of  the 
test  apparatus.  This  is  referred  to  as  system  compliance.  Compliance  was 
found  to  be  a  function  of  the  load  and  can  be  removed  from  the  resulting  data. 
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IV.    EXPERIMENTATION 

The  experimental  setup  consisted  of  a  Material  Tester  (INSTRON  4206) 
and  an  integrated  data  acquisition  system  consisting  of  an  Instron  Control 
Console,  IBM  PC/AT  and  software,  a  detail  description  and  a  list  of  procedures 
is  in  Appendix  B.  The  object  of  the  experiment  was  to  test  a  current  composite 
fiber  type  (used  in  many  Navy  aircraft)  with  a  range  of  several  gauge  lengths 
and  attempt  to  extract  the  underlying  distribution  of  the  fiber.  The  sample 
preparation  is  described  in  Appendix  B.  Each  bundle  sample  was  tested,  the 
data  stored,  and  quick  analysis  conducted  to  ensure  no  improper  testing  method 
were  corrupting  the  data. 

A.      BUNDLE  TESTING 

The  fiber  material  tested  was  made  from  a  Hercules  Magnamite  high 
strength  graphite,  type  AS-4  spool  145  in  3k  bundles  (i.e.  there  are  3,000 
filaments  in  a  bundle).  Several  different  gauge  lengths  were  tested  in  order  to 
analyze  the  effects  of  friction,  slack,  compliance  on  the  fiber  size.  Size  affect 
(Appendix  A)  should  be  able  to  correlate  the  location  parameter  beta  between 
fiber  lengths  if  weakest  link  physics  is  applicable.  In  all,  nineteen  bundle 
samples  where  tested.  The  test  procedure  was  conducted  as  specified  in 
Appendix  B.    The  complete  test  log  is  listed  in  Table  I. 

The  bundle  samples  were  handled  with  extreme  care  so  as  not  to  damage 
the  any  of  the  fibers.  Moving  the  bundles  from  assembly  bench  to  the  storage 
bench  and  then  to  the  testing  machine  presented  the  greatest  hazard.  The 
difficulty   came   when   picking   up   the   sample.      This   was    done   manually   and 
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Table  I. 

TEST  LOG. 

Sample 

#                G.L. 
(cm) 

Oil/ 
Dry 

Comments 

090901 

5.0 

dry 

Good 

090902 

0.5 

dry 

No  failure-slip 

090903 

0.5 

dry 

No  failure-slip 

090904 

0.5 

dry 

No  failure-slip 

090905 

0.5 

dry 

No  failure-slip 

090906 

0.5 

dry 

No  failure-slip 

090907 

5.0 

dry 

Slip 

090908 

5.0 

dry 

Good 

090909 

5.0 

oil 

Good 

090910 

2.5 

dry 

Possible  damage 

100901 

2.5 

dry 

Good 

100902 

2.5 

oil 

Slip 

100903 

2.5 

oil 

Slip 

100904 

50.0 

dry 

Good  -  (friction) 

100905 

50.0 

oil 

Good  -  (friction 
reduced) 

100906 

50.0 

oil 

Possible  damage 

100907 

25.0 

oil 

Good 

100908 

25.0 

dry 

Good 

bending  the  bundles  slightly  was  unavoidable.  The  long  bundles  (500  and  250 
mm)  were  less  vulnerable  to  damage  because  of  the  relatively  large  aspect  ratio 
but  these  were  difficult  to  steady  during  transportation  and  tended  to  bounce 
around. 

Placing  the  sample  in  the  grips  also  presented  operational  difficulties  when 
tightening  the  upper  grip  with  a  socket  wrench.  An  attempt  was  made  to 
manually  counter  torque  the  grip,  but  this  tended  to  be  jerky  at  best.  Also  the 
grip  pressure  was  difficult  to  control  between  samples.  The  alignment  of  the 
bundle  in  the  grips  was  another  concern.  This  was  done  by  sight.  If  not  done 
properly  the  bundle  failure  mechanism  might  be  influenced.  Manually  reducing 
the  experimental  slack  (not  to  be  confused  with  bundle  slack)  was  very  touchy. 
When   toggling  the   grips   at   the  control   console   the   load   read   out  had   to   be 
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acutely  monitored  for  the  bundle  loading  must  be  kept  to  a  minimum  prior  to 
testing. 

During  testing  a  real  time  graphical  display  of  the  load  vs  displacement 
was  used  for  monitoring  the  results,  with  this  experimental  irregularities  in  the 
several  of  the  samples  was  detected.  The  bundles  in  the  very  short  gauge 
lengths  showed  irregular  load  curves.  Keying  into  this  after  several  samples 
helped  to  identify  a  problem  with  the  sample  preparation. 

Friction  was  expected  to  affect  the  long  gauge  length  (250  and  500  mm) 
and  not  affect  the  shorter  gauge  lengths  (25  and  50  mm).  To  test  the  affect 
each  gauge  length  tested  a  pair  was  treated  with  oil  and  a  pair  was  left 
untreated  (dry).  Tests  were  run  on  samples  of  500  mm,  250  mm,  50  mm, 
25  mm  and  5  mm. 

B.      DATA  REDUCTION 

After  tests  were  complete  the  data  files  were  decoded  into  ASCII  and 
stored  for  processing.  The  data  output  from  the  Instron  software  is  a  three 
column  file.  The  first  column  is  the  record  number,  the  second  is  the 
displacement  (inches)  and  the  third  is  the  load  (Ibm).  This  is  then  processed 
through  a  program  to  convert  the  data  to  mm  and  kg  and  remove  compliance. 
There  are  two  data  files  output,  one  file  with  compliance  removed  and  the  other 
converted  raw  data  (compliance  not  removed)  [Ref.  5].  The  output  from  the 
conversion  program  is  then  processed  through  an  analysis  program  listed  in 
Appendix  C.  This  outputs  the  ECDF,  the  bundle  Modulus  E  and  the  Slack 
region  ECDF.  The  ECDF  is  then  processed  through  a  Weibull  Maximum 
Likelihood    Estimator    program    to  define  the  shape    parameter  alpha    and  the 
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location  parameter  beta.    ECDF  is  also  processed  through  a  linearizing  program 
which  outputs  the  Weibull  probability  plot  data  file. 
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V.  RESULTS 

A.      TEST  OBSERVATIONS 

During  testing  several  of  the  specimens  were  invalid  due  to  experimental 
difficulties.  Some  of  the  samples  had  maximum  strains  noticeably  beyond 
expected  limits.  One  bundle  was  observed  to  have  a  much  lower  modulus  than 
what  was  expected  from  equation  C.5.  Dissection  of  the  samples  with  suspect 
loading  curves  revealed  a  lack  of  adhesive  wetting  resulting  in  a  large 
percentage  of  the  fibers  slipping  from  the  tab.  The  problem  appeared  most 
frequently  in  the  small  gauge  length  (GL)  Samples.  It  is  suspected  that  the 
adhesive  which  couples  the  bundle  and  tab  was  not  wetting  the  internal  fibers 
of  the  bundle,  with  a  resulting  loss  of  ability  to  carry  the  higher  loads  of  the 
small  GL.  This  phenomenon  manifested  itself  two  ways.  If  there  was  partial 
slipping  after  initiation  of  the  failure  region,  a  strain  higher  than  the  actual 
maximum  strain  was  observed.  If  total  slipping  occurred  in  a  large  group  of 
fibers  the  modulus  decreased.  After  closer  inspection  of  all  the  samples  it  was 
seen  that  many  of  the  specimens  had  a  small  number  fibers  slipping  in  the  tabs 
but  the  percentage  was  small  enough  not  to  effect  the  results. 

To  avoid  this  problem  the  remaining  tests  were  conducted  in  samples  at  or 
beyond  25  mm.  The  problem  was  less  acute  in  the  long  GL  because  of  the 
lower  loads  these  long  GL's  experience  as  a  consequence  of  size  effect. 

Two  types  of  test  were  carried  out,  with  oil  and  dry.  This  data  was  then 
reduced  by  categories  of  compliance  removed  and  not  removed.  The  compliance 
was  estimated  by  testing  samples  of  zero  GL  [Ref.  5].  The  resulting  curves 
were  fitted  and  reduced  to  an  analytical  equation,  displacement  as  a  function  of 
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load.  This  equation  can  then  be  used  to  subtract  the  compliance  from  the 
empirical  data.  The  zero  GL  samples  used  the  same  coupling  adhesive  and 
possibly  suffered  from  partial  slipping.  This  causes  the  compliance  curve 
coefficients  to  be  over  estimated  and  as  a  result  displacements  removed  for 
compliance  were  too  large.  When  the  compliance  was  removed  from  the 
experimental  data  displacements  in  the  lower  range  showed  negative  values 
indicating  over  estimations. 

Since  the  determination  of  compliance  is  slightly  in  error  it  is  not  used  in 
the  results  discussion,  though  data  reduction  calculations  were  completed  to 
indicate  possible  trends  and  affects.  Table  II  shows  the  data  reduction  with 
removal,  Table  III  is  without  compliance  removal  and  Table  IV,  for  comparison, 
is  the  result  from  single  fiber  tests  [Ref.  6]. 
Table  II.    EXPERIMENTAL  DATA  RESULTS,  COMPLIANCE  REMOVED. 


Compliance  Removed 

Test  # 

Gage 

Modulus 

alpha 

Beta 

dry/oil 

length 

(kg/mn 

0 

(mm) 

(gm/e) 

(mm/mm) 

100901 

25 

109.76 

915 

3.0 

0.0157 

dry 

090901 

50 

64.09 

1068 

3.214 

0.0135 

dry 

100908 

250 

11.0 

917 

3.18 

0.00972 

dry 

100904 

500 

5.433 

905 

3.035 

0.0076 

dry 

090909 

50 

54.7 

912 

4.465 

0.0145 

oil 

100907 

250 

11.18 

932 

3.964 

0.0101 

oil 

100905 

500 

5.755 

959 

3.66 

0.0085 

oil 
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B.      DATA  COMPARISON 

Generally  the  results  indicate  only  small  deviations  in  the  alpha  parameter 
(Table  III)  for  all  GL.  Between  a  oiled  sample  of  the  same  GL  only  the  500 
mm  GL  had  any  noticeable  change  in  the  shape  parameter,  which  had  been 
expected.  The  location  parameter  was  effected  most  by  the  oil  treatment,  beta's 
increased  by  10%  in  both  the  50  and  500  mm  GL.  In  250  mm  GL  the  alpha 
and  the  beta  parameters  seem  to  be  relatively  stable,  meaning  the  oil  treatment 
had  little  effect. 

It  was  expected  that  friction  would  affect  the  longer  GL's  and  have  little 
or  no  affect  on  the  shorter  samples.  The  500  mm  test  did  show  a  considerable 
difference  when  oil  was  applied  (Figure  5). 

The  friction  was  expected  to  become  influential  only  as  the  percentage  of 
failed  fibers  increased,  but  the  dry  bundle  shows  a  decreased  load  early,  which 
can  be  seen  in  the  ECDF  plot  (Figure  6).  In  fact  the  oil  and  dry  sample  data 
merge  as  the  percentage  of  failures  increases.  The  reduced  friction  made  a 
difference  in  the  resulting  Weibull  parameters,  alpha  increased  from  3.085  to 
3.8.  From  single  fiber  testing  alpha  (Table  IV)  for  AS-4  is  approximately  4.11. 
Even  with  oil  the  bundle  is  affected,  as  seen  by  the  sudden  drops  in  the  bundle 
load  (Figure  5)  or  the  sharp  increase  in  the  percentage  of  failures  (Figure  6)  at 
discrete  points. 

Without  friction  the  failure  distributions  are  expected  to  be  smooth. 
Decreasing  the  GL  it  is  anticipated  that  the  effects  of  friction  diminishes.  Tests 
at  250  mm  the  loading  curve  is  seen  to  smooth  out  (Figure  7),  providing  a 
slightly  more  continuous  plot.  Friction  had  less  effect  on  this  GL  though  a 
slight  decrease  is  noticed  when  the  bundle  is  dry. 


19 


Table  III.    DATA  REDUCTION,  COMPLIANCE  NOT  REMOVED. 


Compliance 

Not  Removed 

Test  # 
dry/oil 

Gage 
length 
(mm) 

Modulus 
(kg/mm) 

(gm/e) 

Alpha            Beta 
(mm/mm) 

Beta 
sized  to 
50  mm 

090901 
dry 

100908 
drv 

100904 
dry 

50 

250 

500 

45.6 
9.85 
5.2 

761 
820 
867 

4.75 
4.07 
3.085 

0.01615 

0.0104 

0.00781 

0.01615 

0.0154 

0.0137 

090909 

50 

43.5 

725 

4.89 

0.0179 

0.0179 

on 
100907 

250 

10.6 

833 

3.95 

0.01041 

0.0154 

on 

100905 

oil 

500 

5.5 

917 

3.8 

0.00868 

0.0152 

SUMMARY 

a 

Std 

B              Std. 

E 

Std 

OIL                  4.21        0.59        0.0162      0.0015 
DRY                3.96        0.84        0.0151      0.00125 
Merged           4.09        0.66        0.0157      0.00137 

825 
816 
821 

96 
53 
75 

Table  IV.    RESULTS  FROM  SINGLE  FIBER  TEST  AT  50  mm  GL. 

Fiber  Load  Weibull  Parameters 

Spool  #  Alpha  Beta 

008  4.28  0.0175 

019  3.94  0.0172 

average  4.11  0.0173 

Modulus  -  900  gm/e 
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The  ECDF  plot  (Figure  8)  indicates  friction  had  more  influence  than  was 
revealed  in  the  load  curve.  The  upper  tails  again  merge,  indicating  either  the 
oil  no  longer  decreases  the  friction  or  friction  no  longer  affects  the  remaining 
fibers  in  this  region  of  the  ECDF.  A  slight  change  can  be  seen  between  the 
two  samples  in  the  lower  and  mid  region  but  this  could  be  variations  between 
bundles  and  not  whether  the  specimen  was  treated  with  oil,  but  it  should  not 
be  dismissed  because  it  is  consistent  with  the  other  samples. 

A  50  mm  GL  test  (Figure  9)  showed  results  which  had  been  anticipated  in 
the  longer  lengths.  Little  change  is  observed  as  few  fibers  fail  but  the  affects  of 
friction  become  prominent  as  the  percentage  increases. 

The  ECDF  (Figure  10)  indicates  the  effect  of  oil  more  clearly.  The  entire 
range  of  failure  is  shifted  not  only  at  the  upper  tail  but  the  lower  tail  as  well. 

Determining  how  well  the  underlying  Weibull  distribution  is  extracted  from 
a  bundle  test  and  what  influence  the  oil  treatment  has  can  be  observed  in 
Weibull  plots.  For  the  500  mm  sample  (Figure  11  &  12)  friction  has  definitely 
affected  the  distribution,  but  it  is  also  observed  that  oil  helped  to  reduce  some 
of  the  deviation  especially  in  the  lower  tail. 

As  the  gauge  length  is  reduced  it  becomes  apparent  from  the  Weibull  plots 
(Figures  13,14,15  &  16)  for  the  250  mm  GL  and  the  50  mm  GL  that  friction  is 
not  necessarily  negligible,  for  the  distribution  is  affected.  This  is  clearly 
demonstrated  in  the  50  mm  GL  treated  with  oil  (Figure  15)  which  shows  the 
characteristic  Weibull  distribution  distinctly  as  apposed  to  the  dry  sample. 
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Figure  7.    250  mm  gage  length  load  curve  comparison. 
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Figure  10.    ECDF  for  50  mm  comparison. 
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Figure  15.    50  mm  GL  Wiebull  plot,  oil  treated. 
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Figure  16.    50  mm  GL  Wiebull  plot,  dry. 
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VI.  SUMMARY  AND  RECOMMENDATIONS 

A.  SUMMARY 

The  affect  of  friction  is  evidenced  in  all  gauge  lengths  and  is  not  always 
limited  to  the  upper  tail.  Lubrication  of  the  bundle  reduces  the  effect  of  friction 
and  should  be  used  in  all  bundle  failure  testing. 

The  compliance  removal  was  not  used  in  the  data  reduction  due  to  the 
error  in  the  compliance  curve  caused  by  adhesive  coupling  failure  in  the  tabs. 
The  compliance  ended  up  being  over  estimated.  This  a  side  the  trend  of 
compliance  removal  on  the  data  reduction  (Table  II)  indicated  a  general 
reduction  in  the  parameters  (cc,G).  Estimating  the  compliance  to  be  much  less 
and  extrapolating  this  thought,  it  can  be  assumed  the  value  for  the  parameters 
presently  determined  would  be  slightly  less.  In  a  comparison  of  the  single 
fibers  parameters  with  the  bundle  parameters  (Table  III),  the  bundle 
parameters  were  barely  higher.  Given  the  affect  of  compliance  it  can  be 
concluded  that  once  removed  the  bundle  parameters  will  be  indistinguishable 
from  the  single  fiber  testing  and  in  fact  given  the  quantity  of  the  statistical 
base  much  more  reliable. 

It  was  shown  that  friction  effected  all  the  gauge  lengths  but  the  longer 
gauge  length  (500  mm)  experienced  the  most  disturbance.  The  50  mm  gauge 
lengths  when  treated  with  oil  produced  a  very  smooth  distribution  in  the 
Weibull  plot.  The  dry  50  mm  sample  produced  a  poor  Weibull  plot.  The  50 
mm  gauge  length  produce  good  results  through  out  the  entire  range.  The  500 
mm  gauge  length  produced  good  results  in  the  Weibull  plot  at  the  lower  tail 
(Figure  7).    Since  the  lower  tail  of  the  large  gauge  length  are  statistics  of  much 
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higher  reliability  than  the  50  mm  the  linear  statistics  of  the  500  mm  can  be 
pieced  into  the  50  mm  statistical  plot  there  by  extending  the  range  of  the 
variable  (strain).  The  upper  tail  of  the  50  mm  gauge  length  is  not  completely 
linear  in  the  Weibull  but  this  can  be  removed  and  statistics  from  a  5  mm  test 
can  replace  it.  This  piece  meal  work  on  the  Weibull  plot  can  be  done  because 
of  size  effect.  The  slope  alpha  does  not  change  but  beta  decreases  in  the 
Weibull  as  the  gauge  length  increases.  With  this  procedure  friction  and  other 
contaminates  can  be  removed  and  the  underlying  parameters  accurately 
predicted. 

Through  simulation  a  normal  or  uniform  distribution  of  slack  was 
determined  to  have  very  little  effect  on  the  distribution  if  the  dO  displacement 
is  subtracted. 

B.     RECOMMENDATIONS 

To  ensure  better  results  several  areas  of  experimentation  can  be  improved 
on. 

-  The  handling  of  the  specimens  should  be  reduced,  a  rigid  carrying 
mechanism  which  can  transport  the  bundle  from  assembly  to  storage 
to  testing  without  any  human  handling  or  disturbing  the  fibers  is 
necessary. 

-  The  manual  mechanical  grips  should  be  replaced  with  a  pneumatic  device. 

-  Oil  treatment  should  be  used  on  all  samples. 

-  Improve  the  tab  coupling. 

The  slack  distribution  requires  more  investigation  than  was  addressed  here. 
Simulation  was  a  key  element  which  helped  to  determine  it's  effect  on  the 
distribution.  An  analysis  algorithm  is  provided  in  this  thesis  which  will  help  to 
determine   the   slack   ECDF.      To   further   study   the   question  of  slack   several 
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bundle  samples  (large  gauge  length)  dedicated  solely  for  the  determination  of 
slack  should  be  tested.  A  much  larger  data  base  is  required  to  use  the  methods 
discussed  in  Appendix  C.  By  decreasing  the  cross-head  rate  and  increasing  the 
acquisition  rate  a  large  data  file  can  be  attained.  Once  the  distribution  is  fitted 
to  a  model,  simulation  can  determine  what  affect  it  has  on  the  bundle  strength 
distribution. 
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vn.  CONCLUSIONS 

The  results  achieved  in  this  thesis  conclude  convincingly  that  the 
distribution  of  the  fiber  can  be  retrieved  through  bundle  failure  tests.  The 
parameters  extracted  from  the  testing  compared  favorably  with  the  single  fiber 
testing.  As  a  result,  enormous  characterization  efficiency  is  realized;  in  fact 
with  this  technique,  statistics  to  106  are  realized. 
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APPENDIX  A 
SIZE  EFFECT 

The  determination  of  the  scale  parameter  beta  of  a  gage  length,  which  is  too 
short  to  preclude  practical  experimental  implementation,  can  be  accomplished  by 
size  effect  of  the  weibull  distribution  function. 

Starting  with  the  Weibull  CDF  for  a  fiber  of  defined  gage  length  lx: 

F,(x)  =  1  -  expMx/B,)"1}  (A.l) 

R,(x)  =  1  -  F,(x)  (A.2) 

For  a  second  fiber  of  gage  length  12  with, 

1,  <  12    and    \2f\]  =  m  ;      m  =  integer  >  0 


and  assuming  the  reliability  R,  is  constant. 


R,  =  (RI)(R2XR3)....(RJ  =  R>m  (A.3) 


substituting  A.2  for  R, 


Rt  =  expKx/B.r}"1  (A.4) 
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and  the  total  reliability 


Rt  =  expMx/BJH  (A.5) 


equating  equation  A. 4  and  A.5 


expWx/BJ-}    =    expMx/VT  (A.6) 


taking  the  natural  log  A.6  reduces  to 


(x/Bmr    =  mCx/B,)"1  (A.7) 


again  taking  the  natural  log  of  A.7 


aJn(x/BJ  -  aJnMx/B,)}  =  0  (A.8) 


since  o^,  =  o^  =  a  is  not  =  0 


6,/B2  =  m1/a  (A.  9) 
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APPENDIX  B 
EXPERIMENTAL  PROCEDURES 

A.  EXPERIMENTAL  SET  UP 

The  Test  apparatus  consisted  of  an  INSTRON  Model  4206  material  tester, 
with  a  50.0  kg  load  cell,  connected  to  a  4200  Series  Expanded  Control  Console. 
The  Console  was  connected  to  a  IBM  PC/AT  through  com  port  1  via  a  IEEE 
connection  and  converter.  The  data  acquisition  and  material  tester  control  was 
commanded  by  Instron  Series  DC  Automated  Material  Testing  System  series 
4.01C  software. 

The  controlling  software  has  a  5  page  menu  which  requires  setting  prior  to 
testing.  Many  of  these  settings  are  for  Instron  data  reduction  processing  which 
are  not  used  in  this  thesis  but  required  values  to  operate.  Tho  settings  crucial 
to  testing  and  proper  acquisition  are  cross-head  speed,  data  acquisition  rate  and 
gage  length.  Safety  features  include  a  maximum  load  and  displacement  setting 
to  insure  the  load  cell  and  testing  model  are  not  damaged. 

B.  INSTRON  CALIBRATION  AND  TESTING 

The   Instron  4206   and  control  unit  required   a  warm   up  time  of  one   hour 
prior  to  operation  to  allow  for  stabilization.     The  load  cell  has  two  methods  for 
calibration,  an  electronic  calibration  and  a  mechanical  calibration.    A  mechanical 
calibration  was  performed  for  all  tests. 
1.    Mechanical  Calibration 

Any  inputs  refer  to  the  4200  Series  Expanded  Control  console. 
-    Main  power  on/off/on  -  wait  for  diagnostic 
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-  Press  Load  Balance  /  Enter 

-  Hand  prescribe  a  5.0  kg  weight 

-  Press  Load  Cal 

-  Enter  the  weight  /  Enter 

-  Remove  the  weight 

-  Enter  Load  Balance  /  Enter 

-  Re-hang  the  weight  to  verify  correct  calibration 

-  Repeat  if  required 

2.  Loading  the  Sample 

The    procedures    for   sample    placement   into    the    test    apparatus    are    as 
follows. 

-  The  Instron  specimen  grips  are  separated  enough  to 
allow  the  sample  to  hang  from  the  top  grip  without 
touching  bottom  grip.    The  tab  of  the  specimen  is 
firmly  clasped  with  tweezers  and  gently  lifted  off 
resting  bench  allowing  to  hang  free. 

-  The  samples  top  tab  is  placed  in  the  top  grip, 
careful  to  align  the  bundle  with  the  center,  which 
has  been  measured  and  marked.    After  alignment  the 
grip  is  tightened  making  sure  not  to  twist  or  damage 
the  specimen. 

-  A  load  balance  is  performed  at  the  control  panel  to 
remove  the  weight  of  the  specimen. 

-  The  cross-head  is  then  toggle  down  so  as  to  place 
the  bottom  tab  in  the  bottom  grip  then  tighten. 
Note:  after  tightening  the  grips  the  specimen  develops 
noticeable  slack  and  the  load  cell  measures  a  negative 
load.    The  specimen  fibers  are  slightly  compressed  and 
need  to  be  straightened. 

-  The  cross-head  is  toggle  up  until  a  load  of  0.1  kg 
is  indicated  in  the  load  readout  window  on  the 
control  panel.    The  load  is  then  brought  down  to 
zero  and  the  displacement  reference  zeroized.    The 
procedure  is  repeated  once  to  ensure  the 
displacement  returns  to  at  or  near  zero.    The 
specimen  is  ready  for  testing. 

3.  Testing 

The   test  is   initiated   at  the   IBM   PC/AT   through   the   Instron   software. 
After  the  specimen  test  parameters  are  set  and  loaded,  prior  to  test  initiation, 
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the  IEEE  port  is  enabled.  The  software  then  signals  the  system  is  ready  to 
begin  testing. 

As  the  cross-head  displaces,  load  and  displacement  data  are  sent  to  the 
PC/AT  at  regular  intervals  of  displacement,  as  prescribed  by  the  initial  software 
settings. 

A  real  time  graphical  output  is  available  to  monitor  the  testing.  This 
allowed  for  a  quick  identification  of  bad  tests.  After  the  test  was  complete,  the 
data  file  is  converted  from  Instron  System  code  to  ASCII  for  latter  data 
analysis. 

C.  SAMPLE  PREPARATION 

The  samples  were  made  from  a  Hercules  Magnamite  high  strength  graphite, 
type  AS-4,  spool  145.  The  bundle  has  3000  fibers  with  a  denier  of  .0057446 
gm/in.. 

1.    Procedures 

A  length  of  bundle  is  clamped  at  one  end  on  a  3  meter  aluminum  track 
and  weighted  on  the  other.  Half  Copper  tabs  approximately  2.54  cm  (width)  x 
2.54  cm  (length)  are  pre-positioned  on  the  track  at  a  specified  gage  length. 
Slots  built  into  the  track,  center  the  bundle.  An  adhesive  is  applied  at  the  tab, 
the  bundle  and  the  remaining  half  of  the  copper  tab  are  securely  clamped  on  to 
complete  the  tab.  The  bundle  is  severed  between  samples  on  the  track  and  the 
adhesive  is  allowed  to  cure. 

Once  dry  the  samples  are  manually  removed  from  the  track  and  placed 
on  the  sample  bench. 
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D.  CALCULATION  OF  CROSS-HEAD  SPEED 

The  cross-head  calculation  is  a  function  of  the  number  of  data  points 
required,  data  acquisition  rate,  load  cell  and  displacement  gage  tolerance. 

The  time  to  bundle  failure  was  selected  to  remain  constant.  This  will 
provide  a  constant  strain  rate  on  the  different  gage  lengths  that  were  tested. 

The  relationship  between  the  cross-head,  acquisition  rate  and  recorded  range 
is  given  in  equation  B.l. 

XH  =  (RRXARVDPTS  (B.l) 

where:       XH  =  cross-head  speed 
RR  =  recorded  range 
AR  =  acquisition  rate 
DPTS  =  number  of  data  points  recorded 

The  sampling  interval  is 

deltaD  =  RR  /  DPTS  =  XH  /  AR  (B.2) 

In  order  not  to  sample  beyond  the  abilities  of  the  acquisition  system  the 
following  rule  must  be  followed 

deltaD  >  Recording  Tolerance  (B.3) 

The  Instron  output  tolerance  for  displacement  was  0.00254  mm,  therefore 
deltaD  >  0.00254  mm.  Using  this  as  a  guideline  and  the  constant  testing  time 
chosen  of  5  min..    Table  IB  lists  the  Instron  test  setting  calculations. 
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Table  IB.    TEST  SETTINGS. 


Gage 

length 
(mm) 

Max 
displ. 
(mm) 

Time 
(min) 

Cross- 
head 

(mm/s) 

Acq. 
rate 
(pts/s) 

500 

25 

5 

8.33xl02 

6.6 

250 

12.5 

5 

4.17xl02 

3.3 

50 

2.5 

5 

8.33xl03 

3.3 

25 

1.25 

5 

4.17xl03 

1.67 

44 


APPENDIX  C 

DATA  ANALYSIS 

The  analysis  of  the  bundle  failure  load  vs  displacement  curve  can  be  divided 
into  three  distinct  regions  (see  Figure  1C): 

-  Slack  region. 

-  Linear  region. 

-  Non-linear  region. 

Slack  is  created  during  specimen  production  and  placement  in  the  Material 
testing  machine.     As  the  bundle  is  measured,  cut  and  tabs  are  placed  on  each 


LOAD  CURVE 
REGIONS 


failure 


linear 


slack 


mniTiiTMirili i  ii  "i  i 8iJ "~  -  - 


Figure  1C.    Regions  of  the  load  curve. 

end  of  the  sample  some  of  the  individual  fibers  become  loose  and  effectively 
longer  than  other  fibers.  This  produces  an  uneven  loading  across  the  fibers  as 
the  cross-head  of  the  tensile  testing  machine  strains  the  bundle.     The  ilh  fiber  is 
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not  loaded  until  the  displacement  of  the  cross-head  is  equal  to  the  respective 
slack  of  the  ith  fiber. 

When  the  cross-head  displaces  beyond  the  maximum  slack,  the  loading  curve 
becomes  linear  until  the  first  failure. 

The  amount  of  slack  varies  according  to  an  undefined  distribution.  The 
slack  distorts  the  failure  region  thereby  disturbing  the  underlying  sequence  of 
the  strength  distribution.  The  effect  slack  has  on  the  underlying  shape 
parameter  depends  on  the  slack  distribution  and  amount  of  slack  relative  to  the 
total  strain,  this  is  demonstrated  in  Appendix  D. 

A.    SLACK  DISTRIBUTION 

The  slack  in  a  bundle  is  dependent  on  the  variation  in  the  fiber  lengths. 

L,  =  1  +  delLj  (C.l) 

where:       Lj  =  length  of  fiber  j 
1    =  mean  gage  length 
delLj  =  slack  of  fiber  j 

Since  the  delLj  «  1  the  length  is  approximated  by 

h}  =  1  (C.2) 

The  load  in  the  slack  region  can  be  defined  by  Hooke's  law. 

k-l 

Ps(dA)k  =  I  E  [(dA  -  delLjVd  +  deity]  (C.3) 


k-l 

-  I  E  d/1  ;       k=l,2,....,n 

j=i 

where:       Ps(  )k  =  bundle  load  as  each  fiber  k  is  loaded 
E  =  fiber  modulus  (gm/strain) 
dA  =  apparent  displacement  =  d^  +  delL, 
dV  =  fiber  displacement 
n  =  number  of  fibers 
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'  figure  2C.  Histogram  representation  of  a  PDF  plot 


Defining  the  loading  slack  in  terms  of  the  distribution  function  (fs)  of  Ps  at 
a  displacement  dA  as  shown  in  Figure  2C.    Ps  can  be  write  as 


Ps(d)1  =  {fsCd^Xd^O)  E  U(d-dA1){(d-dA1)/l} 


(C.4) 


Al 


where:  d  =  any  displacement  such  that  d  >  d 

U(d)  =  unit  step  function  (=1  d  >  dA1; 

{=0  d  <  dA1 

Summing  the  load  (eqn.  C.4)  over  the  range  of  dA    and  defining  the  bundle 

modulus  as  a  function  of  the  fiber  modulus 


E  =  nE/1 


Ps(d)  =  E  I  [fs(dAiXdAl+1  -  dAi)l  U(d-dAl) 


(C.5) 
(C.6) 
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differentiating  with  respect  to  d  for  a  continuous  function 


fd, 


6P/5(d)  =  E  fs(d)  d(d)  =  E  F(d,)  (C.7) 


F(dj)  is  the  CDF  of  the  slack.    Assuming  the  empirical  data  approximates 
the  continuous  data. 


P*  =  F 


F  =  P' 


(C.8) 


The   derivative   of  the   empirical   data   can  be   numerically   solved   with   two 
methods,  a  numerical  differentiation  or  a  discretization. 

1.     Numerical  Differentiation 

Even  interval  Forward  and  Central  difference  methods  are  used  due  to 
the  nature  of  the  slack  loading  and  the  recording  system.  As  the  fibers  load 
the  slope  is  discontinuous  as  seen  in  Figure  3C.  Several  numerical  methods 
exist      which      will      produce     varying 


accuracy.  In  trying  to  recover  the 
underlying  distribution  of  the  slack 
the  sampling  rate  and  number  of 
fibers  must  be  considered  (i.e.,  how 
many  data  points  can  be  recorded 
relative  to  the  number  of  fibers). 

If  the  sampling  rate  is  fast 
and  the  number  of  total  data  points  i 
recorded  in  the  slack  region  is  2-3 
times    that    of    the    fiber,  a    3  point 


Slack  Region 


(RepresenUtion) 


J 
I 


•   underlying  loadinj 
D   discrete  data 


figure  3C.  Slack  region  showing 
discrete  data  on  the  underlying  loading 
slope. 
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forward  difference  method  may  yield  accurate  results  due  to  the  linearity 
between  discrete  data.  If  the  relative  sampling  rate  is  slow  the  3  Pt.  forward 
difference  may  skew  the  distribution.  To  prevent  this  a  central  difference 
method  is  safe. 

3  pt.  Central  Difference  method. 

PT  =  (-EH  +  Pui)/2h  (C.9.1) 

5  pt.  Central  Difference  method. 

E,' =  (Ei.2-8P,r8Pi+1-Pi+2)/2h  (C.9.2) 

where:  h  =  $lAi+]  -  <iA, 

2.     Discretization 

The  data  in  the  slack  region  can  be  fit  to  a  curve  and  a  equation 
estimated.  This  equation  can  be  differentiated  analytically.  The  underlying 
data  itself  is  not  smooth  so  such  a  method  is  not  recommended. 

B.    LINEAR  REGION 

The  slope  of  the  linear  region  is  the  bundle  modulus.  To  calculate  this 
slope  a  least  squares  method  is  very  accurate,  due  in  part  to  the  natural 
linearity  of  the  data.  Before  this  method  can  be  employed  the  region  must  be 
specifically  denned.  The  point  where  the  non-linear  slack  region  ends  is 
identified  as  dl  and  the  end  of  the  linear  region  is  identified  as  d2  (the  point  of 
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the  first  failure).     This  is  graphically  displayed  in  Figure  4C.     Between  dl  and 
d2  the  least  squares  method  calculates  modulus  E. 

The  points  dl  and  d2  must  be  graphically  identified  and  manually  entered 
into  the  analysis  program.  To  do  this  a  plot  of  the  load  curve  must  be 
generated  before  the  analysis  can  be  conducted. 


Figure  4C.    The  linear  region. 

1.     Least  Squares  Method 

Due    to    the    noise,    the    discrete    data   i    will    deviate    slightly    from    the 
underlying  slope  E.    From  Hooke's  law  the  expected  relationship  is 


E  =    =  ED, 


where:      D]    =  dArdl 


(CIO) 
(C.10.1) 
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?igure  5C.    Least  Squares  Method 


But    due    to    noise    the    load    P    deviates    by    an    amount    Pe,    this    is 
graphically  represented  in  Figure  5. 


Pe,  =    ED,    -    P, 


(C.ll) 


Defining  the  least  squares  sum  M 


M(E)  =  I  {(ED,  -  PJ/2G2}2 


(C.12) 


where:  a  =  standard  deviation 
Differentiating  M  with  respect  to  E  and  setting  it  equal  zero  defines  the 
minimum  least  squares. 

dM/dE  =  0 
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substituting  in  C.10.1  this  reduces  to 

E  =    IP^-dD/K^-dD2  (C.13) 

2.    Intercept  of  the  Abscissas 

With  the  Modulus  E  the  intercept  of  the  abscissas  can  be  calculated. 

dO    =  dl  -  P(dl)/E  (C.14) 

C.    FAILURE  REGION 

This  region  is  the  key  to  retrieving  the  strength  distribution.  Using  the 
slope  of  the  bundle  the  failure  statistics  can  be  transferred  to  a  ECDF  space. 
From  the  ECDF  the  parameters  of  the  distribution  can  be  estimated. 

The  ECDF  is  a  plot  of  the  percentage  of  failures  versus  a  variable. 
Typically  this  variable  is  the  strength,  this  would  be  difficult  to  translate  due  to 
the  nature  of  the  bundle  test.  The  displacement  is  constantly  incremented, 
where  as  the  load  is  a  bundle  load  and  can  vary.  Because  of  this  the  ECDF  is 
plotted  versus  the  displacement.  The  displacement  is  proportional  to  the 
strength  as  defined  by  Hooke's  law  and  once  the  fiber  modulus  is  identified  the 
displacement  can  be  translated  to  load. 

When  a  fiber  fails,  the  bundle  modulus  E  changes  by  one  fiber  strength. 
From  Hooke's  law  it  is  expected  that  the  failures  would  occur  along  the 
modulus  Ej.  Ej  being  the  modulus  of  the  unfailed  bundle  which  ideally  is  a 
multiple  of  the  single  fiber  modulus  proportional  to  the  number  of  remaining 
fibers.  Unfortunately  this  is  not  the  case,  as  the  bundle  loads,  friction  and 
slack  change  the  distribution  causing  the  actual  loading  path  to  deviate  from  E} 
and  follow  along  E,. 
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To  show  this  in  more  detail  consider  first  the  ideal  underlying  loading  (no 
slack,  friction  or  noise)  of  a  bundle  under  constant  strain.  At  a  break  in  one 
fiber  j,  the  load  will  drop  until  intercepting  the  slope  of  E,.,  shown  in  Figure  6C. 
This  new  slope  emanates  radially  out  from  the  intercept  of  the  abscissas,  (note: 
without  slack  the  intercept  should  be  the  origin).  If  the  data  acquisition  rate 
where  high  enough  this  failure  could  be  clearly  identified.  Now  add  slack,  this 
extends      the      failure 


displacement  of  some 
fibers  and  the  intercept 
of  the  abscissas  by  E, 
after  the  initial  failure 
to  some  point  between 
the  origin  and  the  dO 
point.  The  intercept 
will  shift  a  relatively 
small  amount  but  will 
digress  toward  the 
origin.         If    friction     is 


D  Failure 
j   |  fibers 

remaining 
—  underlying 

loading 

•    discrete 

data 


^igure  6C.    Failure  loading  path. 


present  the  broken  fibers  will  cling  to  the  unbroken  fibers  and  cause  additional 
strain.  Add  noise  to  the  acquisition  system  and  the  discrete  data  fluctuates. 
These  factors  combine  to  force  the  actual  recorded  data  to  deviate  from  the 
underlying  bundle  loading  path. 

A  bundle  with  n  number  of  fibers  will  have  n  failures.  If  the  data 
acquisition  system  can  not  record  this  many  points  in  the  failure  displacement 
range  than  each    statistic    can  not  be  recorded.      With  slack,  friction    and  noise 
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even   with   a   fast   enough    data   acquisition   system   the   failures   might   not   be 
distinguishable. 

The  percentage  of  failures  in  the  bundle  test  can  be  defined  by  the  ratio  of 
the  modulus  of  no  failure  to  the  modulus  of  failures.  The  reliability  is  defined 
by 

E,  =  E/E  =  EOU/tE^-dO)]  (C.15) 

The  percent  failure  is  1  minus  reliability,  the  ECDF  is  defined  by 


EOL)  =  l  -  E, 


(C.16) 


This  relationship  is  graphically  depicted  in  Figure  6C.    The  ECDF  begins  at 


Figure  7C.    Transfer  from  loading  curve  to  ECDF  space. 
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d2,  so  the  manual  selection  of  this  point  must  be  done  with  care. 

Assuming  all  conditions  ideal  the  ECDF  would  fall  out  as 

F  =  j/n  (C.17) 

But  because  conditions  are  not  ideal  this  does  not  occur.  What  is  important  is 
the  shape  of  the  ECDF  and  whether  it  is  close  to  the  underlying  CDF.  The 
ability  to  distinguish  each  failure  is  unimportant.  Another  point  to  note  is  the 
lower  tail  needs  to  be  clearly  defined.  Due  to  the  relatively  low  failure  rate  in 
the  lower  tail  and  a  constant  recording  rate  this  area  of  the  curve  has  a  high 
number  of  discrete  data  points  per  failure,  so  in  essence  the  testing  method 
takes  care  of  itself. 

1.    Upper  Bound 

The  data  is  recorded  on  even  intervals  of  displacement  not  just  at  the 
failure  points,  this  is  graphically  described  in  Figure  6.  The  discrete  data  most 
closely  resembling  the  failure  point  is  on  the  top  edge  of  the  loading  curve.  This 
is  referred  to  as  the  Upper  bound.  When  this  is  translated  into  the  ECDF 
space  the  Upper  bound  is  on  the  far 
edge  of  the  plot  as  seen  in  Figure  8C. 
To  retrieve  the  underlying  CDF  from 
the  ECDF  the  Upper  bound  (Fu) 
must  be  extracted. 

Attempting  to  extract  the 
Upper  bound  from  a  bundle  of  n 
fibers  the  most  failure  points  that 
could  be  optimistically  retrieved  is  n. 
With  a  large  n  this  is  not  possible 
due     to     the     acquisition     rate,     no     Figure  8C.    Close  view  of  the  ECDF 
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matter,    the    shape    and    location    is    important    and    does    not    require    a    large 
number  of  points. 

Using  the  empirical  ranking  as  a  reference  point 

F  =  k/m  ;      k  =  1,2, ,m  (C.18) 

where  m  is  the  number  of  points  desired  in  the  ECDF.     F  is  used  to  define 
search  bands  of  data  where  the  upper  and  lower  band  is  defined  by     1/n.    The 
value  of  F,  which  has  the  largest  displacement  (o^)  is  assigned  to  Fuk.     When 
no  discrete  data  is  in  the  search  band  F,  Fuk  is  assigned  the  value  of  F  at  a 
displacement  dA  k.,. 

2.     Maximum  Likelihood  Estimator. 

With  the  Upper  bound  extracted  this  data  can  be  used  to  estimate  the 
parameters  of  the  distribution.  The  model  selected  is  a  Two  Parameter  Weibull 
Distribution.  The  estimator  chosen  is  the  Maximum  Likelihood  Estimator 
(MLE)  [Ref.  #:p.  103]  developed  by  R.  A.  Fisher.  This  method  finds  the 
maximum  of  the  likelihood  of  the  sample  in  the  range  of  the  variable  o!u.  and  is 
a  function  of  distribution  parameters.  For  the  Weibull  model  the  method  solves 
for  the  two  parameters  alpha  and  beta.  The  parameter  beta  is  solved  explicitly 
[Ref.  4]  in  terms  of  alpha.  Alpha  can  not  be  solved  explicitly  and  an  iteration 
must  be  used. 

13  =  [  1/m  I  du,2  ]m       sum  i=l  to  m  (C.19) 

a  =  {  I  du,2  ln(du,)/(I  du,2)  -  I  lnCdu,)}-1  (C.20) 
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An  Atkinson  method  of  iteration  was  chosen,  initially  alpha  is  assumed.     From 
this  successive  alpha's  are  calculated  using  equation  C.20. 

ct'0  =  a0  -  (a0  -  a,)2/(a0  +  2a,  +  aj  (C.21) 

where:       a0  =  initial  start  value 

a,  =  f0(aj 
Oa  =  fa(a,) 

This  equation  is  iterated  for  a'0  until    a'0  =  a0. 

1.    Linearized  ECDF 

The  data  is  linearized  for  a  Weibull  probability  plot  to  visually  check 
the  conformity  to  this  model.  To  linearize  the  ECDF  start  with  the  Weibull 
equation 

F  =  1  -  exp(-P/fl)a 
subtracting  1,  multiplying  by  -1  and  taking  the  natural  log. 

ln(l-F)  =  (-P/B)a 
again  taking  the  natural  log 

ln[-ln(l-F)l  =  cdn(P/13)  (C.22) 

C.22  is  a  linear  line  with  a  slope  of  alpha. 

D.   DATA  OUTPUT 

As  the  analysis  program  calculates  certain  data  file  are  output,  these  files 
are: 

-  EXPER.OUT,  the  converted  experimental  output  displacement,  strain 
and  load. 

-  SLACK.OUT,  the  slack  region. 

-  EXMS.OUT,  the  load  curve  minus  slack  region. 
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UBFECDF.OUT,  the  Upper  bound  ECDF. 
FECDF.OUT,  the  ECDF. 
SECDF.OUT,  slack  ECDF. 
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APPENDIX  D 

A.      SIMULATION 

The  simulation  of  a  bundle  test  has  two  purposes.  First  by  analyzing  in 
depth  the  mechanisms  which  govern  the  failure  of  the  fibers  in  a  bundle,  a 
better  understanding  of  how  to  extract  the  characteristic  properties  is  attained. 
The  properties  desired  are  the  distribution  parameters. 

Secondly  a  simulation  gives  us  the  ability  to  know  what  the  underlying 
parameters  are  and  thereby  proof  the  analysis  data  reduction  algorithm.  Also 
the  variables  which  influence  the  shape  of  the  failure  distribution  can  be  varied 
in  simulation  and  the  degree  of  induced  affect  observed. 

The  bundle  loading  curve  can  be  divided  into  three  regions; 

-  Slack 

-  Linear 

-  Failure 

Initially   the   simulation   requires   input   of  the   variable   estimates,   this   is 
done  through  an  input  file. 
1.      Initial  Inputs 

The  underlying  parameters  define  the  strength  distribution  for  a  gage 
length.  The  distribution  chosen  is  the  Weibull  model.  The  parameters  are 
alpha  and  beta,  the  number  of  fibers  per  bundle  is  n.  The  fiber  modulus  (E)  is 
in  gm/mm2,  from  Hooke's  law 

o  =  E  e  (D.l) 


where:  o  =  stress 

e  =  strain 
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Substituting  the  definition  of  stress  and  strain  for  fiber  i 


P/A  =  E  6/L,  (D.2) 

where:  P  =  load 

A  =  area  =  constant 

L  =  length  of  fiber 

8  =  change  in  length  L 

rearranging  D.2  and  letting  Ef  =  EA  (gm/e) 

P,  =  Ef  5/L,  (D.3) 

2.       Slack  and  Strength  Distributions 

The  slack  is  caused  by  the  variation  in  fiber  length  L.  If  L  were 
constant  no  slack  would  be  evident.  The  slack  is  expected  to  have  a  definable 
distribution,  but  this  is  unknown  and  could  vary  from  bundle  to  bundle 
therefore  several  models  are  simulated.    Two  models  selected  where 

-  Two  parameter  Weibull 

-  Uniform. 

Although  physical  consideration  favors  the  slack  model  to  be  normal, 
the  Weibull  model  was  chosen  for  it's  analytical  convenience.  The  Weibull  is 
very  flexible  and  can  characterize  a  wide  range  of  shapes.  For  an  alpha  of  3.5 
the  distribution  is  approximately  normal,  other  than  3.5  the  shape  of 
distribution  is  skewed  negative  or  positive.  Beta  identifies  the  central  tendency 
of  the  distribution,  its  relatable  to  the  mean.  From  equation  2.1  the  Weibull 
CDF  is  defined  as 

Fw(x;cc,fl)  =  1  -  exp{-(x/B)a 
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rearranging  in   terms   of  the   variable   x,,   which   can   represent  slack  (delL,)  or 
strength  P,. 

x,  =  exp{  [ln(-ln(l-FJ  +  oln(B)]  /  a)  (D.4) 

0  <  Fw  <  1 

The  uniform  model  is  a  simple  two  variable  distribution  with  a  lower 
(u,)  and  upper  (u2)  range.  Since  the  slack  is  initiated  at  zero  the  lower  range  is 
zero. 

Fu(delL;u„u2)  =  (delL  -  u,)/(u2  -  u,)  (D.5) 

0  <  Fu  <  1 
solving  for  delL,  from  D.5  and  setting  Uj  =  0 

delL,  =  u2  Fu  (D.6) 

3.       Data  Generation  and  Continuity 

A  uniform  random  number  generator  creates  an  array  of  numbers 
between  0  and  1  and  stores  it  until  recalled. 

Recalling  equation  D.3  and  substituting  in  C.l  and  solving  for  the 
displacement  (6,). 

5,  =  P,  (1  +  delL*)  /  Ef  (D.7) 

The  bundle  is  put  under  a  constant  strain  rate,  so  the  discrete 
variable  will  be  displacement.     As  the  bundle  is  strained  the  filaments  will  fail 
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in  order  of  the  underlying  strength  or  strain  (neglecting  slack  and  friction).  For 
a  constant  filament  length  the  strain  is  proportional  to  the  displacement  and 
the  fiber  will  fail  at  some  underlying  displacement.  If  slack  is  present  in  the 
bundle  each  fiber  will  not  fail  at  the  underlying  displacement  because  the 
displacement  is  no  longer  proportional  to  strain  but  a  function  of  both  slack  and 
stain. 

The  displacement  observed  by  the  Instron  material  testing  machine  is 
not  5.    The  machine  records  an  apparent  displacement  (dA) 

dA ,  =  5,  +  delL,  (D.8) 

This   will   affect   the   fiber  failure   order,   every  other  variable   is  then 
keyed  to  the  order  (n,)  of  dA  s. 
4.       Slack  Loading 

As  the  bundle  is  strained,  the  fibers  are  loaded  by  a  progressive 
summation  of  each  fiber  in  order  of  their  slack  distribution.  The  loading 
function  in  the  slack  region  can  be  defined  as  a  function  of  the  displacement  (5.) 
each  fiber  experiences  after  slack  is  released. 


j-i 
P.,  =    ZE,  5,/L,  ;j=  l,2,....,n  (D.9) 

i  =  l 


j-l 

P1J=    LEr  [(delL0j  -  delLo.Vd  +  delL0  ,)1  (D.10) 

i  =  l 


where:  delLo ,  =  order  slack 
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This  relationship  can  be  represented  graphically  in  Figure   ID.     The  last  load 
calculated  from  equation  D.10  will  occur  at  delLon. 


Ps3 
PS1 


delL      delL  n       delL, 

1  2  I 


Figure  ID.    Slack  region  load  simulation. 

4.      Failure  loading 

The  order  of  each  fiber  failure  is  determined  by  the  ordered  strength 
distribution  (assuming  no  slack).  With  slack  the  order  of  failure  will  deviate 
slightly.  The  failure  load  is  a  function  of  the  modulus  of  the  bundle  and  the 
displacement  of  each  fiber.  The  bundle  failure  load  (Pu)  can  be  defined  from 
equation  D.2  as 


PUJ  =  IEf  S/L, 


(D.ll) 


substituting  in  C.l  and  D.8 


PUJ  =  IE  [(dAj  -  dell^/d  +  delLj 


(D.12) 


>=j 
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5.       Discrete  Data 

The  material  testing  machine  measures  and  transmits  data  on  even 
intervals  of  displacement.  To  simulate  this  the  defined  loading  points 
(P„delLo)  and  (Pu,dA)  will  be  set  limits  for  a  linear  interpolation.  The  data  is 
interpolated  on  intervals  of  deldA  which  is  set  by  the  cross  head  speed. 
Defining  the  discrete  apparent  displacement 

^M«    flU,  +  deldA  (D.13) 

The  interpolation  algorithm  to  define  the  discrete  data  from 
continuous  is  the  same  for  both  slack  and  linear  region. 

The  discrete  load  P  in  the  slack  region  is  a  function  of  P.,  the  ordered 
slack  (delL0)  and  the  apparent  discrete  displacement 

E  =  KP.i.i  "  P.jMdelL^  -  delLoPld*.  (D.14) 

The  discrete  data  is  defined  this  way  until  the  first  underlying 
displacement  dA ,. 

The  slope  of  the  linear  region  intercepts  the  abscissas  at  a  point  dO. 
The  slope  of  the  linear  region  is  the  bundle  modulus  E    (see  Figure  2C). 

E  =  (P.,  -  P.n)/(dA1  -  delLJ  (D.15) 

The  slope  can  be  solve  using  E 

dO  =  delL,  -  Ptn/E  (D.16) 
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figure  2C.    Linear  calculations. 

Subsequent  discrete  data  points  can  be  created  by  a  regressive 
process.  At  each  fiber  loss  a  new  slope  emanates  radial  outward  from  dO. 
Between  breaks,  to  define  the  discrete  data,  the  slope  is  defined  by  the  failure 
load  Pu,  the  displacement  dA  and  dO 


£  =  [Pu/(dAj  -  dO)]  (iL.-dO) 
where:         <IA ,      dAj  ;  j=l,2,....,n 


(D.17) 


6.       Noise 

Noise  was  simulated  simply  by  setting  a  tolerance  level  and  using  the 
random  number  generator  to  produce  a  set  of  random  numbers  between  -0.5 
and  0.5.    This  data  was  added  to  each  discrete  load  P. 
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The  input  file  to  the  simulation  has  units  of  kg,mm  and  the 
experimental  output  from  the  Instron  Software  is  in  lbm,in.  A  conversion  of  the 
data  is  required,  this  is  accomplished  through  a  separate  program 

B.      SIMULATION  INTERPRETATION 

The  simulation  was  written  in  fortran  code,  listed  in  Appendix  F.  Multiple 
simulations  were  run  varying  the  parameters  slack  and  noise  to  see  the  effect. 

Simulation  data  was  analyzed  via  the  algorithm  as  represented  in 
Appendix  C.  The  output  is  the  bundle  modulus  E,  the  shape  parameter  of 
distribution  Alpha,  and  the  location  parameter  Beta  (mm/mm). 

The  input  parameters  were  alpha  =  4.0,  beta  =  0.16  kg  and  a  50.0  mm 
Gage  length  (GL).  The  noise  was  varied  from  zero  to  about  0.1%  fluctuation  of 
the  maximum  load.  The  maximum  slack  input  was  3%  of  the  gage  length  and 
a  normal  slack  distribution.  The  seed  used  for  the  random  number  generator 
was  held  constant  so  the  underlying  strength  distribution  was  constant.  The 
seed  for  the  slack  region  and  noise  was  not  held  constant.  The  results  of 
multiple  simulation  test  is  in  Table  IC.  A  value  of  4.134  for  the  shape 
parameter  was  retrieved  consistently  (same  seed),  as  expect  no  deviation  for  the 
same  data.  The  increase  in  slack  drove  the  parameter  down.  The  noise  did  not 
have  any  effect  on  the  parameters.  Figure  3C  show  the  loading  curve  with 
slack  (minus  dO)  and  no  slack. 

No  noticeable  change  in  the  parameter  output  nor  in  the  ECDF  plot 
(Figure  4C). 
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Table  IC.        SIMULATION    OBSERVATIONS    AND 
PROOF  TESTING  RESULTS. 


Maximum  Slack 

alpha 

beta 

o% 

GL 

(kg) 

(13  = 

a/2) 

0 

4.134 

0.016 

5 

4.010 

0.016 

10 

4.099 
4.095 
4.091 
4.089 
4.088 
4.075 
4.076 

0.016 

The  numerical  derivative  of  the  slack  region  is  very  sensitive  to  the  any 
noise.  This  can  be  seen  in  the  plot  of  the  slack  ECDF  Figure  5C.  In  order  to 
retrieve  the  slack  distribution  the  recording  system  must  have  a  low  tolerance 
load  cell. 
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Figure  3C.    Simulation  load  curve,  slack;  no  slack. 
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Figure  4C.    The  simulation  ECDF  for  slack  and  no  slack. 
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Figure  5C.    Simulation  of  the  slack  ECDF.    Weibull  approximation  of  Normal. 
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APPENDIX  E 

I.  DATA  ANALYSIS  CODE 
A.   ANALYSIS  PROGRAM 

1.  Data  Input 

The  program  requires  a  formatted  data  file  called  EXPER.OUT.  The 
data  is  read  in  and  prompts  the  user  for  several  variables.  These  are  listed 
below. 

-  What  is  the  gage  length? 

-  How  many  fibers  in  the  bundle? 

-  Which  Numerical  differentiation  is  desired? 
3,5,7  pt.  forward  cliff,  or  3,5  pt.  central  diff. 

-  The  end  of  the  slack  region  and  the  beginning  of  the 
failure  region  in  displacement  (mm)  (this  requires  a  graphical 
estimate). 

-  How  many  points  in  the  input  data  file? 

-  The  acquisition  tolerance  of  the  displacement?  (this  is  for 
any  filtering  of  data  points  in  the  event  the  sampling  rate 
is  too  high  and  repeated  data  is  evidenced. 

-  The  modulus  is  calculated  and  displayed  and  then  ask  if  this 
value  needs  to  be  change.    This  is  in  the  event  the  modulus 
calculation  in  correct  a  manual  entry  can  be  made.  (Be  Careful 
not  to  corrupt  the  true  data). 

2.  Data  Output 

The  data  out  put  is  to  the  screen  and  to  a  file.  The  following  is  a 
list  of  the  file  outputs. 

-  EXPER.OUT,  the  input  file  is  output  with  an  added  column  of 
strain 

-  SLACK.OUT,  the  slack  region  (load  vs  displacement) 
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-  EXMS.OUT,  the  load  curve  minus  slack  region 

-  UBFECDF.OUT,  the  Upper  bound  ECDF 

-  FECDF.OUT,  the  ECDF 

-  SECDF.OUT,  the  slack  region  ECDF 

B.  MAXIMUM  LIKELIHOOD  ESTIMATOR 

The  MLE  program  reads  one  of  the  ECDF  data  files,  users  choice,  the  first 
line  in  the  input  file  is  a  utility  line  defining  the  file  size,  the  number  of 
filaments  and  the  bundle  modulus  E,  the  remaining  lines  are  data. 

1.  Data  Input 

The  input  other  than  the  read  file  is  prompted  from  the  screen.     The 
following  is  a  list  of  requested  information. 

-  How  many  points  for  the  Expected  Ranking  (The  MLE  requires  a 
probability  ranking) 

-  First  approximation  for  alpha 

2.  Data  Output 

The  parameter  estimations  alpha  and  beta  (mm/mm)  are  output  to  the 
screen. 

C.  LINEARIZED  PLOT 

This  program  takes  the  ECDF  and  does  a  Weibull  linearizing  (equation 
C.20). 
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II.  FORTRAN  CODE 
A.      ANALYSIS  PROGRAM 

PROGRAM  ANALY 

*  BUNDLE  FAILURE  ANALYSIS 

* 

*  Outline  of  program: 

*  Read  in  the  data  points 

*  NLA  -  Number  In  Array 

************************************************************** 

*  Find  the 

*  -linear  region 

*  -slope  of  the  linear  region  =  EB  (E  BAR) 

*  -intercept  point  of  the  linear  region  with 

*  the  abscissas  =  DO    (d  zero) 

*  -Maximum  slack  or  minimum  linear  displacement  -  DPI 

*  -Maximum  displacement  of  the  linear  region  -  DP2 

* 

*  Input  the  range  of  the  linear  region  by  manually  retrieving 

*  the  denning  points  DPI  and  DP2. 

* 

*  OUTPUT:       SLACK.OUT     -  slack  region. 

*  SECDF.OUT      -  slack  ECDF. 

*  EXMS.OUT       -  experimental  data  minus  slack  region. 

*  FECDF.OUT      -  failure  region  ECDF. 

*  UBFECDF.OUT  -  upper  bound  of  the  failure  region  ECDF. 

* 

INTEGER  NLA,NISA 
PARAMETER(NLA=6000>NISA=1000) 

REAL  D(NLA),P(NIA),PP(NIA),H,LP(NIA),LD(NIA),FS(NIA),F(NIA) 
+,FUB(NIA),DF(NIA),DFUB(NIA),PF(NIA),PFUB(NIA),TEMP(NIA),GL 
INTEGER  COUNT,NIDPl,NIDP2,NI(NIA),SLCT 

CHARACTER*^  Ql,Q2 

* 

********************************************** * ****************** 

*  Input  analysis  control  parameters. 

***************************************************************** 

PRINT  *,  'What  is  the  gage  length?  (mm)' 
READ  *,  GL 

PRINT  *,  '  How  many  filements?' 

READ  *,  N 

PRINT  *,  '  Which  method  do  you  require  for  reduction  of 

PRINT  *,  '  the  slack  region  to  the  ecdf  space?' 
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PRINT  *,  '  1)  3  PT  forward  diff.  method.' 

PRINT  *,  '  2)  5  PT  forward  diff.  method.' 

PRINT  *,  '  3)  7  PT  forward  diff.  method.' 

PRINT  *,  '  4)  3  PT  central  diff.  method.' 

PRINT  *,  '  5)  5  PT  central  diff.  method.' 

READ  *,  SLCT 

* 

***************************************************************** 

*  Read  the  range  of  the  linear  displacement. 

***************************************************************** 

* 

PRINT  *,  'INPUT  DPI  (the  maximum  slack)' 

READ  *,  DPI 

PRINT  *,  'INPUT  DP2  (the  maximum  linear  displacement)' 

READ  *,  DP2 

******************************************************************** 

*  Find  the  how  many  points  are  in  the  data  file. 
* 

******************************************************************** 

PRINT*,  TIow  many  points  in  the  data  file?' 

READ  *,  COUNT 

******************************************************************** 

*  Find  out  what  the  tolerance  of  the  machine  is. 

*  (displacement) 

******************************************************************** 

PRINT  *,  'What  is  the  tolerance  of  the  displacement  output' 
PRINT  *,  'data?  This  is  the  limit  between  data  displacement' 
PRINT  *,  'Input  in  mm.' 
READ  *,  TOLD 

PRINT  *,  TIow  many  points  in  the  empirical  rank?' 
READ  *,  NER 
144  PRINT  *,  T)o  you  want  the  load  in  Kg  (K)  or  Newtons  (N)?' 

READ(*,303)  Q2 
IF(Q2  .EQ.  IT  .OR.  Q2  .EQ.  TSP)  THEN 

PRINT  *,  '  ' 
ELSE 

PRINT  *,  Try  again  in  UPPER  CASE  KORN1 
GOTO  144 

ENDIF 

* 

******************************************************************** 

*  Read  in  the  experimental  results. 

* 

******************************************************************** 

* 

OPEN(UNIT=ll,FILE='EXPER.OUT,STATUS='OLD') 

* 

DO  10  I=l,COUNT 
READUU10)  NI(I),TEMP(I),P(I) 
10       CONTINUE 
REWINDUl) 
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CONVER  =  1.0 

IF(Q2  .EQ.  W  THEN 

CONVER  =  9.8062 

ENDIF 

DO  44  I=l,COUNT 

P(I)  =  P(I)*CONVER 

WRITE(11,1001)NI(I),TEMP(I),P(I),TEMP(I)/GL 

44       CONTINUE 

******************************************************************** 

*  Filter  out  any  points  with  difference  in  displacement  interval 

*  less  than  the  machine  tolerance  TOLD. 

* 
******************************************************************** 

J=l 

D(l)  =  TEMP(l) 

DO  11  I=2,COUNT 

DIFF  =  TEMPd)-TEMP(I-l) 

IF(DIFF  .LT.  TOLD)  THEN 

PRINT  *,  'SKIPPED  POINT  ',1 
ELSE 

J=J+1 

D(J)  =  TEMP(I) 

P(J)  =  P(D 
ENDIF 
11       CONTINUE 

COUNT  =  J 

* 

*  FIND  COUNT  (NI)  OF  DPI  AND  DP2 

******** 

* 

1=1 
20  IF(D(I)  .LT.  DPI)  THEN 

1=1+1 

GOTO  20 
ELSE 

NIDP1  =  I 

DPI  =  LXD 
ENDIF 
30  IF(D(I)  .LT.  DP2)  THEN 

1=1+1 

GOTO  30 
ELSE 

NIDP2  =  1-1 

DP2  =  LXI-1) 

ENDIF 

* 

******************************************************************** 

*  OUTPUT  THE  SLACK  REGION. 

* 

******************************************************************** 

* 
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0PEN(UNIT=13,FILE='SLACK.0Ur,STATUS=1JNKN0WN') 
DO  500  I  =  1,NIDP1 
WRITE(13,H0)  I,EKD,P(D 
500  CONTINUE 

CLOSEU3) 

******************************************************************** 

*  LEAST  SQUARES  TO  DETERMINE  THE   SLOPE   BETWEEN  DPI   AND 

DP2. 

* 

******************************************************************** 

* 

TOP  =  0.0 
BOT  =  0.0 
DO  40  I=(NIDP1+1),NIDP2 

LP(I)  =  P(I)  -  P(NIDPl) 

LLXI)  =  LXI)  -  DPI 

TOP  =  TOP  +  LP(I)*LD(I) 

BOT  =  BOT  +  LLXI)**2 
40  CONTINUE 

EB  =  TOP/BOT 

********* 

*  Find  the  intercept  point  do 
*******  ** 

* 

DO  =  DPI  -  P(NIDP1)/EB 
PRINT  *,  'Do  =  \D0 
PRINT  *,  '  ' 

PRINT  *,  'The  Modulus  calculated  from  the  load  curve' 
PRINT  *,  'is  \EB 
PRINT  *,  '  ' 

PRINT  *,  "Would  you  like  to  change  this  (Y/N)?' 
READ(*,303)  Ql 

IF(Q1  .EQ.  T  .OR.  Ql  .EQ.  'y')  THEN 
PRINT  *,  'Input  the  new  modulus' 
READ  *,  EB 
ENDIF 

303         FORMAT!  Al) 

* 

* 
******************************************************************** 

*  OUTPUT  THE  EXPERIMENTAL  DATA  MINUS  SLACK. 

* 
******************************************************************** 

OPEN(UNIT=13,FILE='EXMS.OUr,STATUS='UNKNOWN') 
J=l 

WRITE(13,H0)  J,0.0,0.0 
DO  550  I  =  NIDP1.COUNT 
J=J+1 

WRITE(13,H0)  J,(D(I)-D0),P(I) 
550  CONTINUE 
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******************************* ************************************* 

*  NUMERICAL  DIFFERENTIATION  OF  THE  SLACK  REGION 

*  (using  a  forward  and  central  difference  methods). 

* 
******************************************************************** 

* 

IFCSLCT  .EQ.  1)  GOTO  103 
IF(SLCT  .EQ.  2)  GOTO  105 
IF(SLCT  .EQ.  3)  GOTO  107 
IFCSLCT  .EQ.  4)  GOTO  130 
IF(SLCT  .EQ.  5)  GOTO  150 

********************************** 

*  3  pt  forward  diff.  method. 

********************************** 

* 

103  DO  3  I=1,NIDP1 

*  if  out  of  the  slack  region  leave  the  loop 

IF  (D(I)  .GT.  DPI)  GOTO  51 
H  =  ABS(D(I+1)-D(I)) 

*  3  point  forward  difference  numerical  differentiation. 

*  Generate  Fs(da)  CDF 

PP(I)  =  (-3*P(I)+4*P(I+l)-P(I+2))/(2*H) 
FS(I)  =  PP(I)/EB 
3         CONTINUE 

GOTO  51 

********************************** 

*  5  pt  forward  difference  method. 

********************************** 

105  DO  5  I=1,NIDP1 

*  if  out  of  the  slack  region  leave  the  loop 

IF  (D(I)  .GT.  DPI)  GOTO  51 
H  =  ABS(D(I+1)-D(I)) 

*  5  point  forward  difference  numerical  differentiation. 

*  Generate  Fs(da)  CDF 

PP(I)  =  C-25*P(I)+48*P(I+l)-36*P(I+2)+16*P(I+3)- 
C3*P(I+4))/(12*H) 
FS(I)  =  PP(I)/EB 

5         CONTINUE 

********************************** 

*  7  pt  forward  difference  method 

********************************** 

107  DO  7  I=1,NIDP1 

*  if  out  of  the  slack  region  leave  the  loop 

IF  (D(I)  .GT.  DPI)  GOTO  51 
H  =  ABS(D(I+1)-D(I)) 

*  7  point  forward  difference  numerical  differentiation. 

*  Generate  Fs(da)  CDF 

PP(I)  =  C-147*P(I)+360*P(I+l)-450*P(I+2)+ 
+400*P(I+3)-225*P(I+4)+72*P(I+5)-10*P(I+6))/(60*H) 
FS(I)  =  PP(I)/EB 
7         CONTINUE 
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********************************** 

*  3  pt  central  difference  method 

********************************** 

* 

130  1=1 

H  =  ABS(EKI+1)-D(I)) 
PP(I)  =  (-3*P(I)+4*P(I+l)-P(I+3))/(2*H) 
FS(I)  =  PP(I)/EB 
DO  33  I=2,NIDP1 

*  if  out  of  the  slack  region  leave  the  loop 

IF  (EKD  .GT.  DPI)  GOTO  51 
H  =  ABS(D(I+1)-D(I)) 
PP(I)  =  (-P(I-1)+P(I+1))/(2*H) 

*  Generate  Fs(da)  CDF 

FS(I)  =  PP(I)/EB 

33       CONTINUE 

********************************** 

*  5  pt  central  difference  method 

********************************** 

* 

150  1=1 

*  Use  a  3  pt  forward  difference  method  for  first  couple  pts. 

H  =  ABS(D(I+1)-D(I)) 
PP(I)  =  (-3*P(I)+4*P(I+l)-P(I+3))/(2*H) 
FS(I)  =  PP(I)/EB 
1=2 

H  =  ABS(EKI+1)-D(I)) 
PP(I)  =  (-3*P(I)+4*P(I+l)-P(I+3))/(2*H) 
FS(I)  =  PP(I)/EB 
DO  55  I=3,NIDP1 

*  if  out  of  the  slack  region  leave  the  loop 

IF  (EKI)  .GT.  DPI)  GOTO  51 

H  =  ABS(D(I+1)-D(I)) 
PP(I)  =  (P(I-2)-8*P(M)+8*P(I+l)-P(I+2))/(12*H) 

*  Generate  Fs(da)  CDF 

FS(I)  =  PP(I)/EB 
55       CONTINUE 

51       CONTINUE 

******************************************************************** 

*  FAILURE  REGION 

******************************************************************** 

J=0 

DO  60  I=NIDP2,COUNT 
J=J+1 

F(J)  =  1  -  P(IV(EB*(D(I)-D0)) 
DF(J)  =  EKD-DO 
60       CONTINUE 

NUMF  =  J 

******************************************************************** 

*  EXTRACT  THE  UPPER  BOUND. 
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* 

************************************* ******************************* 

DEL  =  1.0/NER 

K=NUMF 

I=NER 

FUBQ+1)  =  F(K) 

DFUBG+1)  =  DF(K) 

K=K-1 

1=1-1 
FL  =    1.0*I/NER  -  DEI72 
FU  =    1.0*I/NER  +  DEL/2 
PRINT  *,TL  F(K)  FU' 

PRINT  *,  FL,'  <  ',F(KV  <  \FU 

90  IF((K  XT.  1)  .OR.  (I  XT.  0))  GOTO  92 

91  IF((F(K)  XT.  FU)  .AND.  (F(K)  .GE.  FL))  THEN 

FUB(I+1)  =  F(K) 
DFUB(I+1)  =  DF(K) 

1=1-1 
FL  =  1.0*I/NER  -  DEIV2 
FU  =  1.0*I/NER  +  DEL/2 

K=K-1 

IF(K  XT.  1)  GOTO  92 

GOTO  91 
ELSEIF(F(K)  XT.  FL)  THEN 

*  check  to  see  how  many  it  skips 

******** 

Jl  =  1 

93  IF(F(K)  XT.  (1.0*(I-J1)/NER  -  DEL/2))  THEN 

J1=J1+1 

GOTO  93 

ENDIF 

DO  94  I1=I,(I-J1+D,-1 

IF(I1+1  XT.  1)  GOTO  92 

FUB(I1+1)  =  1.0*I1/NER 

DFUB(I1+1)  =  DF(K) 
94  CONTINUE 

1=11 

FUBQ+1  )=F(K) 

DFUB(I+1)=DF(K) 

1=1-1 
FL  =  1.0*I/NER  -  DEL/2 
FU  =  1.0*I/NER  +  DEL/2 

K=K-1 

GOTO  90 
ELSE 
K=K-1 

IF(K  XT.  1)  GOTO  92 
GOTO  90 
ENDIF 

92  CONTINUE 

PRINT  *,  "OUT  OF  LOOP  WITH  K  AND  I  ' 
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PRINT  \  'K=  \K,'    1=  ',1 

NUMEF  =  1+2 

***************** ********************************************* 

*  WRITE  THE  DATA 

* 
************************************************************** 

OPEN  (UNIT=ll,FILE=rUBFECDF.OUr,STATUS=rUNKNOWN') 
MSUM  =  NER+1-NUMEF 
WRITE(  11,1010)  MSUM,N,EB 

DO  250  I=NUMEF,NER 
WRITE(11,110)  I,DFUB(I),FUB(I) 
250  CONTINUE 

CLOSE(ll) 

OPEN(UNIT=ll,FILE=TECDF.OUr,STATUS=TJNKNOWN') 
WRITE(  11,1010)  NUMF,N,EB 
DO  70  I=1,NUMF 
WRITE(11,110)  I,DF(I),F(I) 
70       CONTINUE 
CLOSE(ll) 

OPEN(UNIT=12,FILE='SECDF.OUT,STATUS='UNKNOWN') 
DO  80  I=1,NIDP1 
WRITE(12,110)  I,EKI),FS(I) 

80       CONTINUE 

******************************************************************** 

*  FORMATS 

* 
******************************************************************** 

100  FORMAT(lX,I5) 

110  F0RMAT(1X,I5,2X,F8.4,2X,F8.4) 

1001  FORMAT(1X,I5,2X,F8.4,2X,F8.4,2X,E10.4) 

1000  FORMAT(1X,I5,2X,E16.10,2X,E16.10) 

1010    FORMAT(  1X,I5,2X,I5,2X,F8.4) 

************************************************************ 

END 

B.      MAXIMUM  LIKELIHOOD  ESTIMATOR 

PROGRAM  MLE 

INTEGER  NIA 

PARAMETER  (NIA  =  6000) 

REAL  FE(NIA),DFE(NIA),FUB(NIA),DFUB(NIA),GL 

CHARACTER*1  Ql 

* 

PRINT  *,  'Do  you  want  the  upper  bound  or  complete  data' 
PRINT  *,  'set,  suggest  the  compete  data  set  when  the' 
PRINT  *,  'DATA  is  less  the  the  number  of  elements' 

110         PRINT  *,  TYPE  "U"  for  upper  and  "C"  otherwise' 

303         FORMATYA1) 

READ  (*,303)  Ql 

IF(Q1  .EQ.  'U)  THEN 

OPEN(UNIT=  1 1  ,FILE='UBFECDF.OUr,STATUS='OLD') 
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ELSEIF  (Ql  .EQ.  'O  THEN 

0PEN(UNIT=11,FILE='FECDF.0UT\STATUS='0LD') 

ELSE 

PRINT  *,  TRY  TYPING  CAPITAL  LETTERS' 

GOTO  110 

ENDIF 

READ(11,1010)  M,N,EB 
1010       F0RMAT(1X,I5,2X,I5,2X,F8.4) 

DO  10  I=1,M 

READ(11,1020)  J,DFUB(I),FUB(I) 
10  CONTINUE 

1020       F0RMAT(1X,I5,2X,F8.4,2X,F8.4) 

************************************************************** 

*  CALCULATE  ALPHA  AND  BETA 

PRINT  *,  rWhat  is  the  GAGE  LENGTH?' 

READ  *,  GL 

* 

*  REDUCE  THE  DATA  TO  ABOUT  50  POINTS 
888         PRINT  *,  How  many  points  for  the  MLE  ?' 

READ  *,  MLE 
PRINT  *,  MLE 
DEL  =  1.0/(MLE+1.0) 
J  =  1 
I  =1 

*  find  the  value  of  F  above  and  below  the  expected  rank  and  average 
400         FE(I)  =  1.0*I/(MLE+1.0) 

IF(I  .GT.  MLE)  GOTO  450 
410         IF((FUB(J)  .LE.  FE(D)  .AND.  (FUB(J+1)  .GT.  FE(I)))  THEN 

DFE(I)=  -(DFUB(J+1)-DFUB(J))*(FUB(J+1)-FE(I))/(FUB(J+1)-FUB(J)) 
+    +  DFUB(J+1) 
1=1+1 

GOTO  400 

ELSEIF((J+1)  .LE.  (M))  THEN 
J=J+1 
GOTO  410 
ENDIF 
450         PRINT  *,  I,J,FUB(J),FUB(J+1),FE(I) 
IF(FUB(J+1)  .LE.  0.0)  THEN 
PRINT  *,  'MLE  to  large  try  again' 
GOTO  888 

ENDIF 

************************************************************** 

*  estimate  the  initial  alpha 

PRINT  *,  'Estimate  the  intial  alpha' 

READ  *,     ALFA0 

************************************************************** 

*  calculate 

ALFA  =  ALFA0 
KK=0 
320  DO  310  IT  =  1,2 

SUM1=0.0 
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SUM2=0.0 
SUM3=0.0 
*  Calculate  the  variables  of  alpha 
DO  300  1=1, MLE 

IF(DFEd)  .GT.  0.0)  THEN 

SUM1  =  SUM1  +  ((DFE(I))**ALFA  *  LOG<DFE(I))) 
SUM2  =  SUM2  +  ((DFE(I))**ALFA) 
SUM3  =  SUM3  +  LOCXDFE(I)) 

ENDIF 
300  CONTINUE 

IF(IT  .EQ.  1)  THEN 

ALFA1=1.0/(SUM1/SUM2  -  (SUM3/MLE)) 

ALFA  =  ALFA1 
ELSEIFdT  .EQ.  2)  THEN 

ALFA2=1.0/(SUM1/SUM2  -  (SUM3/MLE)) 
ENDIF 
310  CONTINUE 

ALFA0P  =  ALFA0  -((ALFA0  -  ALFA1)**2)/(ALFA0  -  2*ALFA1 
+  +  ALFA2) 

TEST  =  ABS(ALFAOP-ALFAO) 
IF(TEST  .GT.  .0001)  THEN 

IF(KK  .GT.  1000)  GOTO  341 
KK  =  KK  +  1 
ALFAO  =  ALFAOP 
PRINT  *,  'ALFAO  =',ALFA0 
GOTO  320 
ELSE 

ALPHA  =  ALFAOP 
ENDIF 
SUM  =  0.0 
DO  340  I  =  1,MLE 
SUM  =  SUM  +  (DFE(I))**ALPHA 

340  CONTINUE 

341  BETA  =  (SUM/(MLE))**(1/ALPHA) 

PRINT  *,  'ALPHA  =  '.ALPHA,'      BETA  =  '.BETA/GL,'  mm/mm' 
END 

C.      WEffiULL  LINEAR  CDF 

PROGRAM  LIN 
INTEGER  NLA,AREA 
PARAMETER  (NIA  =  6000) 

REAL  FE(NIA),DFE(NIA),FUB(NIA),DFUB(NLA),FS(NIA),D(NIA) 

* 

OPEN(UNIT=ll,FILE=UBFECDF.OUT\STATUS='OLD') 

READUU010)  M,N,EB 
1010       FORMATdX,I5,2X,I5,2X,F8.4) 

DO  10  1=1  ,M 

READU  1,1020)  J,DFUB(I),FUB(I) 
10  CONTINUE 

1020       F0RMAT(1X,I5,2X,F8.4,2X,F8.4) 
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*  REDUCE   THE   DATA   TO  ABOUT   50   POINTS,   USING  THE   EXPECTED 
RANK 

PRINT  *,  HOW  MANY  POINTS  TO  FILTER  USING  THE  EXPECTED 
RANK 

READ  *,  NUM 
PRINT  *,  NUM 
DEL  =  1.0/(NUM+1.0) 
J  =  1 
I  =1 

*  find  the  value  of  F  above  and  below  the  expected  rank  and  average 
400         FECI)  =  1.0*I/(NUM+1.0) 

IF(I  .GT.  NUM)  GOTO  450 
410         IF((FUB(J)  .LE.  FE(D)  .AND.  (FUB(J+1)  .GT.  FE(I)))  THEN 

DFE(I)=  -(DFUB(J+1)-DFUB(J))*(FUB(J+1)-FE(I))/(FUB(J+1)-FUB(J)) 
+    +  DFUB(J+1) 
1=1+1 

GOTO  400 

ELSEIF((J+1)  XT.  (M+1))  THEN 
J=J+1 
GOTO  410 
ENDIF 

450         PRINT  *,  I,J 

******************************************************************** 

*  CALCULATE  ALPHA  AND  BETA 

* 

*  FIND  F*=LOG(-LOG<l-F))  VS  LOG{DELTA/BETA) 

DO  20  I=1,M 

IF(FUBd)  .LE.  0.0)  FUB(I)  =  0.0001 
FS(I)  =  LOG<-LOGU-FUB(I))) 
IF(DFUB(I)  .LE.  0.0)  DFUB(I)  =  0.0001 
EXI)  =  LOGQFUBQ)) 
20  CONTINUE 

DO  70  I=1,M 
510       IF((FS(I)  .LE.  0.63212)  .AND.  (FS(I+1)  .GT.  0.63212))  THEN 
BETA  =  -(D(I+1)-D(I))*(FS(I+1)-0.63212)/(FS(I+1)-FS(I)) 
+    +  EXI+1) 
GOTO  71 
ENDIF 

70  CONTINUE 

71  CONTINUE 
BETA  =  EXP(BETA) 
PRINT  *,  '  BETA  =  '.BETA 


DO  80  I=1,M 

IF(DFUB(I)  .LE.  0.0)  DFUB(I)  =  0.0001 
EXI)  =  LOG<DFUB(I)/BETA) 
80  CONTINUE 

BETA  =  EB*BETA/N 

******************************************************************** 

OPEN(UNIT=  1 1  ,FILE='FSTAR.  OUT  .STATUS^UNKNOWN') 
DO  60  I=1,M 
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WRITEd  1,1020)  I,EKI),FS(I) 
60  CONTESTUE 

END 
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APPENDIX  F 


SIMULATION  FORTRAN  CODE 


PROGRAM  SIMUL 


SIMUL  IS  A  STATISTICAL  SIMULATION  ROUTINE  USED 

FOR  RELIABILITY 

ANALYSIS  OF  COMPOSITE  FILEMENTS. 


UTILITIES  FILES  REQUIRED  ARE: 

SEED  -  file  containing  the  integer  seed  for  random 

number  generation 
INPUT  -  file  containing  the  input  parameters  * 


OUTPUT  FILES  INCLUDE: 

EXPER.OUT  -  the  experimental  simulation  of  discrete 

data. 
SEED  -  outputs  a  new  seed. 


* 
* 
* 
* 
* 
* 
* 
* 
* 
* 
* 
* 
* 
* 
* 
* 
* 
* 
* 
* 
* 
* 
* 
* 
* 
* 

* 
* 
* 
* 
* 
* 
* 
*__ 

*  INPUT  PARAMETERS  ARE: 

*  1)  THE  LENGTH  OF  THE  TEST  BUNDLE  (MGL). 

*  2)  THE  NUMBER  OF  FILEMENTS  PER  BUNDLE  (N). 

*  3)  THE  NUMBER  OF  SEGMENTS  PER  FILEMENT  (M). 

*  4)  THE  SEGMENT  LENGTH  IS  DEFINED  BY  MGL=M*LS 

* 

*  THIS  PROGRAM  CREATES: 

* 


INPUT  FILE  FOR  BUNDLE2.FOR  PROGRAM       :line  1 

XXXXXXXXXXXXXXXXXXX1234567S 

:line  2 

E    (kg) 

.90                           : 

line  3 

ROWS  =    M 

100 

:line  4 

COLS  =    N 

100 

:line  5 

MGL     (mm) 

250.0 

:line  6 

ALPHA 

4.0 

:line  7 

BETA 

16.0 

:line  8 

MODEL 

NORMAL 

:line  9 

SKEW/ALPHA 

3.5 

:line  10 

MAX  SLACK  (% 

)         10.0 

:line  11    %  OF  MAX 
DISPLACEMENT 

SAMPLE  RATE 

20.0 

:line  12    pts/s 

CROSS-HEAD  SPEED         1.8 

:line  13    mm/min 

NOISE    +/-  kg 

.001 

:line  14    tolerance 

TOLD  (mm) 

.00254 

:line  15    tolerance  of  displ 
output 
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*  1)  THE  RANDOM  NUMBERS 

*  2)  CONVERTS  TO  WIEBULL  BUNDLE 

*  3)  FIND  THE  WEAKEST  IN  THE  ENTIRE  BUNDLE 

*  FOR  THE  LOWER  BOUND. 

*  4)  FIND  THE  WEAKEST  IN  EACH  FILEMENT 

*  5)  ORDER  THE  WEAKEST  TO  STRONGEST  FILEMENTS 

*  6)  CHANGE  THE  STRENGTHS  TO  BREAKING  LOADS  AND 

*  DISPLACEMENTS 

*************************************************************** 

*  Begin  Program: 

* 

*************************************************************** 

* 

INTEGER  ROWS,COLS,MAX,LIM,LIMl,N,PTS,COUNT,LIM2 

PARAMETER(LIM=6000,LIM1=1000,LIM2=1001) 

REAL  BUN(LIM1),RAND1(LIM1),RAND2(LIM1),DISPL(LIM1), 
+SLACK(LIM1),P(LIM1),ALPHA,BETA,MGL,BRKLDM> 
+DISPLM,MAXX,MAXY,MEAN,PU(LIM1 ), 
+PL<LIM1),SUM,MAXSLK,TEMP(LIM2),PS(LIM2),D(LIM), 
+DISPLA(LIM1),DEL,PMS,L0AD(LIM),NSE 

DIMENSION  N0(LIM1),NB(LIM1),NDB(LIM1) 

DOUBLE  PRECISION  SEED 

CHARACTERS  MODEL*7,SONS,Q,Q2,Q3,Q4 

COMMON  MAX 

COMMON  /SML/M,N 

* 

******************************************************************** 
* 

*  FILE  11  =  PARAM  DATA  FOR  INPUT 

*  FILE  12  =  SEED  DATA  FOR  SEEDING 

* 

******************************************************************** 
* 

OPEN(UNIT=  1 1  ,FILE='INPUr  ,STATUS='OLD') 

OPEN(UNIT=12,FILE='SEED\STATUS='OLD') 

* 

*  READ  IN  THE  INPUT  DATA 

* 
************************************************************** 

* 

READ(11,1530) 

READ(11,1500)  E 

PRINT  *,  'MODULUS^,  E 

REALX1 1,1520)  M 

PRINT  *,TM  =  \M 

ROWS=M 

REALX1U520)  N 

COLS=N 

PRINT  *,  *N  =  \N 

REALX11,1510)  MGL 

PRINT*,  'Mean  Gage  Length  =  \MGL 
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MAX=M*N 

REATX1U510)  ALPHA 
PRINT  *,  'ALPHA  =  '.ALPHA 
REALXH,1510)  BETA 
PRINT  *,  'BETA  =\BETA 
* 

*  *  *  *  * 

*  Determine  which  model  to  use  for  slack. 


REALX11,1515)  MODEL 

PRINT  *,  TVlODEL^JVfODEL,'  APPROXIMATION7 
IF  (MODEL  .EQ.  'NORMAL')  THEN 
REALX11,1500)  STDDEV 

PRINT  *,  'ALPHA  =  \STDDEV 

REALX11,1510)  MEAN 
PRINT  *,  *Maximum  slack  as  percent  of  expected' 
PRINT  *,  'maximum  displacement  (3%  gage  length). =',ME AN 
MEAN  =  (MEAN/100)*(0.03)*MGL/2 
PRINT  *,  TVlean=',MEAN,'mm' 
PRINT  *,  'Maximum  slack  =',MEAN*2,'mm' 
ELSE 

REALX11,1500)  RANGE1 

REALXH.1500)  RANGE2 
RANGE  1=(RANGE1/100)*(.03)*MGL 
RANGE2=(RANGE2/100)*(.03)*MGL 

PRINT  *,  TIANGE  1  =',RANGE1,'  mm' 

PRINT  *,  TIANGE  2  =',RANGE2,'  mm' 
ENDIF 


*  Read  in  the  sampling  rate  and  cross-head  speed. 
****** 

******************************************************************* 

*  Ensure  the  displacement  interval  governed  by  the  cross-head  and 

*  sampling  rate  is  greater  than  the  tolerance  level  of  the  output 

*  displacement  data. 
* 

* 

READ(11,1510)  SR 
PRINT  *,  'Sampling  rate  =\SR 
READUU510)  XSP 

PRINT  *,  'Cross-head  speed  =',XSP 

* 

****** 

*  How  much  noise. 
****** 

* 

READUU510)  NSE 
PRINT  *,  Tlange  of  noise  =\NSE 
REALXH,1510)  TOLD 
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PRINT  *,  Tolerance  of  output  =',TOLD 

********* 

*  Calculate  the  number  of  total  points  recorded  in  the  * 

*  experiment.  * 

********* 

* 

444      RR  =  MGL*(0.03) 

PTS  =  RR*SR/(XSP/60) 

DELTAD  =  RR/PTS 

IF(TOLD  .GT.  DELTAD)  THEN 

PRINT  *,  The  displacement  interval  is  less  than  the  output' 

PRINT  *,  'tolerance,  please  input  a  new  cross-head  and/or' 

PRINT  *,  'sampling  rate.' 

PRINT  *,  'Current  Cross-head  =  ',XSP 

PRINT  *,  'Input  new  XSP' 

READ  *,  XSP 

PRINT  *,  'Current  sampling  rate  =  ',SR 

PRINT  *,  'Input  new  sampling  rate' 

READ  *,  SR 

GOTO  444 

ENDIF 

* 

*  Determine  the  type  of  simulation. 

****** 

* 

PRINT  *,  'SLACK  OR  NO  SLACK  ?  (Y/N)' 

READ  (*,1600)  SONS 

PRINT  *,  TSfOISE  OR  NO  NOISE  ?  (Y/N)' 

READ  (*,1600)  Q2 

PRINT  *,  'Change  the  SEED  ?' 

READ  (*,1600)  Q 

PRINT  *,  'COMPLAINCE  (Y/N)?' 

READ  (*,1600)  Q3 

IF((Q3  .EQ.  V)  .OR.  (Q3  .EQ.  'y'))  THEN 

PRINT  *,  WHAT  IS  THE  COMPLIANCE  NOISE  LEVEL?' 

READ  *,  CTOL 

ENDIF 

PRINT  *,  T)o  you  want  English  units  (lbs,inches)  on  output? 
+(Y/N)' 
PRINT  *,  'If  N  the  output  is  in  SI  units  (kg,mm)' 

READ  (*,1600)  Q4 

* 

CLOSE(ll) 
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******************************************************************** 

*  BEGIN  THE  ALGORITHM  ** 

******************************************************************** 

** 

*  Read  the  seed,    (seed  for  the  uniform  random  number  generator). 

*  and  generate  the  random  numbers 
********* 

* 

READ(12,1000)  SEED 

REWIMX12) 

DO  10  I=1,N 

CALL  RAND(SEED,M,BUN) 

CALL  SMALL(ROWS,BUN,RANDl(I)) 

10       CONTINUE 

* 

********* 

*  Check  if  slack  is  requested  in  simulation. 
********* 

* 

IF  (SONS  .EQ.  Y)  THEN 
CALL  RAND(SEED,N,RAND2) 

ENDIF 

* 

********* 


*  Initialize  the  order  arrays. 

**  *  *  * 


* 

DO  25  1=1  ,N 
SLACK(I)=0.0 
NO(I)=I 
NB(I)=I 
NDB(I)=I 

25  CONTINUE 

* 

********* 

*  Check  if  seed  is  to  be  changed,  if  so  write  the  new  seed. 

********* 


IF  ((Q  .EQ.  Y)  .OR.  (Q  .EQ.  y»  THEN 
WRITE(  12,1000)  SEED 
REWIND  (12) 
CLOSEU2) 

ENDIF 

* 

******************************************************************** 

*  Generate  the  load  and  slack  distributions. 

* 

******************************************************************** 
* 
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CALL  WEIBUUALPHA,BETA,N,RAND1,P) 

IF  (SONS  .EQ.  V  .OR.  SONS  .EQ.  V)  THEN 
IF  (MODEL  .EQ.  TsTORMAL')  THEN 

CALL  NORMAL(STDDEV,MEANrN,RAND2,SLACK,MAXSLK) 
ELSE 

CALL  UNIFM(RANGEl,RANGE2rN^lAND2,SLACK,MAXSLK) 
ENDIF 
SLACK(1)=0.0 

ENDIF 

* 

*************************************************************+*+**++ 

*  Find  the  displacement  for  each  failure  load  and  slack. 

* 

******************************************************************** 

* 

DO  20  1=1, N 

DISPUI)  =  (P(I)/E)*(SLACK(I)  +  MGL) 
20         CONTINUE 

K=N 

****** 

*  Order  the  "DISPL"  array,  and  carry  "No" 
****** 

* 

CALL  SHELL(N,DISPL,NO) 

* 

****** 

*  reorder  the  remaining  arrays,  the  others  are 

*  keyed  to  "DISPL" 

^p  ^p  ▼  ^  f*  ^p 

* 

CALL  SWITCH(N,NO,RANDl) 
CALL  SWITCH(N,NO,RAND2) 
CALL  SWITCH(N,NO,SLACK) 

CALL  SWITCH(N,NO,P) 

* 

******************************************************************** 

*  Calculate  the  apparent  displacement. 

******************************************************************** 

* 

DO  40  I=1,N 

DISPLA(I)  =  DISPUI)  +  SLACK(I) 

40  CONTINUE 

* 

******* 

*  Reorder  the  arrays  keyed  to  "DISPLA" 
******* 

* 

CALL  SHELL(N,DISPLA,NB) 
CALL  SWITCH(N,NB,RAND1) 
CALL  SWITCH(NJSTB,RAND2) 
CALL  SWITCH(N,NB,SLACK) 
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CALL  SWITCH(N,NB,DISPL) 
CALL  SWITCH(N,NB,P) 
CALL  SWCH2(N,NB,N0) 

******* 

*  Define  the  ordered  slack  as  "TEMP" 
******* 

* 

DO  45  I=1,N 
TEMP(I)=SLACK(I) 
45       CONTINUE 

IF((SONS  .EQ.  Y)  .OR.  (SONS  .EQ.  y»  THEN 
CALL  SHELL(N,TEMP,NDB) 

ENDIF 

* 
******************************************************************** 

*  * 

*  Calculate  the  maximum  slack  and  intercept  point.  * 

*  Calculate  the  upper  and  lower  Pi  bundle  load  at  each  displa.        * 

*  * 

*  * 

*  Find  the  loads  at  the  slack  points. 
******* 

* 

PS(1)=0.0 
DO  70  J=2,N 
CONST=0.0 
DO  75  I=1,(J-1) 

CONST=CONST+E*(TEMP(J)-TEMP(I))/(MGL+TEMP(I)) 
75         CONTINUE 
PS(J)=CONST 

70       CONTINUE 

******************************************************************** 

*  ADD  COMPLIANCE  TO  SLACK  * 

*  * 

*  DC1  =  mm  * 

******************************************************************** 

IF((Q3  .EQ.  T)  .OR.  (Q3  .EQ.  y»  THEN 

DC1  =  .0356 

Al=    0.0512 

A2=  -0.0217 

POWER  =  0.1975 

A3  =  0.0064 

A4  =  0.0242 

DO  120  I  =  1,N 

IF((TEMP(I)  .LT.  DC1)  .AND.  (PS(I)  .GE.  CTOL))  THEN 

TEMP(I)  =  TEMP(I)  +  Al*(PS(I)**POWER)  +  A2 

ELSEIF(PSd)  .GE.  .01)  THEN 
TEMP(I)  =  TEMP(I)  +  A3*PS(I)  +  A4 

ENDIF 
120     CONTINUE 
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ENDIF 

******************************************************************** 

*  FAILURE  POINTS. 

******************************************************************** 

* 

DO  80  J=1,N 

CONST=0.0 
DO  85  I=J,N 

CONST=CONST+E*(DISPLA(J)-SLACK(I))/(MGL+SLACK(I)) 
85  CONTINUE 

PU(J)=CONST 

80  CONTINUE 

******************************************************************** 

*  ADD  COMPLIANCE  TO  FAILURE  REGION 

*  * 

*  DC1  =  mm  * 

******************************************************************** 

IF((Q3  .EQ.  V)  .OR.  (Q3  .EQ.  'y'))  THEN 

DO  121  I  =  1,N 

IF((DISPLA(I)  .LT.  DC1)  .AND.  (PU(I)  .GE.  CTOL))  THEN 

DISPLA(I)  =  DISPLA(I)  +  Al*(PU(I)**POWER)  +  A2 

ELSEIF(PUd)  .GE.  .01)  THEN 
DISPLA(I)  =  DISPLA(I)  +  A3*PU(I)  +  A4 

ENDIF 
121      CONTINUE 

ENDIF 

* 

******************************************************************* 

*  ADD  DISCRETE  DATA  * 

*  TO  THE  * 

*  SLACK  REGION.  * 

*  * 

******************************************************************* 
* 

J=l 

ADD  =  0.0 
CONST=0.0 
LXJ)=DELTAD 
DO  90  I=1,(N-1) 

IF  (EKJ)  .LT.  TEMP(I+1))  THEN 
SLOPE  =  (PS(I+1)-PS(I))/(TEMP(I+1)-TEMP(I)) 
91         LOALXJ)=(D(J)-TEMP(I))*SLOPE  +  ADD 
CONST=LXJ) 
J  =  J  +  1 

LXJ)=CONST+DELTAD 
IF  (EKJ)  -LT.  TEMP(I+1))  GOTO  91 
ENDIF 
ADD  =  ADD  +  PS(I+1)-PS(I) 
90       CONTINUE 
JBIGL  =  J-l 
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******************************************************************* 

*  ADD  DISCRETE  DATA  * 

*  TO  THE  * 

*  LINEAR  REGION  * 

*  * 

******************************************************************* 

* 

CONST=0.0 

SLOPE=(PU(D-PS(N)V(DISPLA(l)-TEMP(N)) 
100         LOAD(J)=SLOPE*(D(J)-TEMP(N))  +  ADD 
CONST=  D(J) 
J=J+1 

D(J)=CONST  +  DELTAD 
IF  (EKJ)  -LT.  DISPLA(l))  GOTO  100 
JENDL  =  J-l 

****** 

*  Find  the  intercept  point. 

****** 

* 

DINT=  DISPLA(l)  -  PU(l)/SLOPE 

* 

******************************************************************* 

*  ADD  DISCRETE  DATA  * 

*  TO  THE  * 

*  FAILURE  REGION.  * 

*  * 

******************************************************************* 

DO  110  I=2,N 

IF  (LXJ)  .LT.  DISPLA(D)  THEN 
SLOPE  =  PU(I)/(DISPLA(I)-DINT) 
115  LOAD(J)=SLOPE*(LXJ)-DINT) 

CONST  =  D(J) 
J=J+1 

EKJ)  =  CONST  +  DELTAD 
IF  (D(J)  .LE.  DISPLA(D)  GOTO  115 
ENDIF 
110         CONTINUE 
LOAD(J)=0.0 
700       PRINT  *,  'CONTINOUS  DATA  COMPLETE  ' 
COUNT  =  J 

PRINT  *,  'COUNT=\J 

* 

******* 

*  Find  the  largest  value  for  plotting 
******* 

* 

CALL  LARGE(K,PU,BRKLDM) 

* 

****** 

*  Where  brkldm  and  sldm  are  the  maximum  for  that  vector. 
****** 
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DISPLM=DISPLA(N) 

* 

****************************************************************** 

*    CREATE  THE  NOISE 

* 

****************************************************************** 


IF  ((Q2  .EQ.  T)  .OR.  (Q2  .EQ.  y»  THEN 
CALL  NOISE(SEED,COUNT,LOAD,NSE) 
ENDIF 

******************************************************************** 

*  PRINT  THE  DATA  * 

*  * 
******************************************************************** 

OPEN(UMT=13,FILE='EXPER.OUr>STATUS=,UNKNOWN') 

IF((Q4  .EQ.  V)  .OR.  (Q4  .EQ.  '/))  THEN 

CONST1=0.03937 

CONST2=2.2046 

ELSE 

CONST1  =  1.0 

CONST2  =  1.0 

ENDIF 

* 

DO  105  1=1, COUNT 

WRITE(  13,1 1 1 )  I,(D(I)*CONSTl ),(LOAD(I)*CONST2) 
105       CONTINUE 

WRITE(13,116)  ALPHA,BETA,N 
WRITE(13,111)  COUNT,DISPLM,BRKLDM 
CLOSE(13) 

*  FORMATS  * 

*  * 

******************************************************************** 

111        FORMAT!  1X,I5,2X,F8.4,2X,F8.4) 

116       FORMAT  (1X,TARAMETERSALPHA=',F4.1,2X,'BETA=', 

+F4.U1X 

+,'N=',I4y,lX,,N,,5X,T)ELTA',10X,'LOAD') 
1000      FORMAT(lX,F15.1) 
1500      FORMAT(20X,E16.10) 
1510      FORMAT(20X,F8.4) 
1515      FORMAT(20X,A6) 
1520      FORMATX20XJ5) 
1530      FORMATUX/) 

1600      FORMATXAl) 

* 

STOP 

END 
*. 

SUBROUTINE  NOISE(SEED,COUNT,LOAD,NSE) 
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INTEGER  COUNT 
DOUBLE  PRECISION  SEED 
REAL  LOAD(COUNT),NSE,RNDM 
DO  10  I=l,COUNT 
CALL  RANDN(SEED,RNDM) 
LOALXD  =  LOAEKD  +  (RNDM*NSE  -  NSE/2) 
10       CONTINUE 
END 


*. 


SUBROUTINE  SWITCH(N,ORDER,ARRAY) 

INTEGER  N,ORDER(N),LIMl 
PARAMETER(LIM=5000) 
REAL  ARRAY(N),TEMP(LIM) 
DO  1  I=1,N 

TEMP(I)=ARRAY(I) 

CONTINUE 

DO  2  1=1  ,N 

K=ORDER(I) 

ARRAY(I)=TEMP(K) 
CONTINUE 
END 


SUBROUTINE  SWCH2(N,ORDER,ARRAY) 
*. 

INTEGER  N.LIM1 
PARAMETER(LIM1=5000) 
INTEGER  ARRAY(N),TEMP(LIMl),ORDER(N) 
DO  1  I=1,N 
TEMP(I)=ARRAY(I) 

1  CONTINUE 

* 

DO  2  1=1  ,N 
ARRAY(I)=TEMP(ORDER(I)) 

2  CONTINUE 
END 
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SUBROUTINE   LARGE(MAX,VECTOR,MAXVAL) 
* 

INTEGER  MAX 

REAL  VECTOR(MAX),MAXVAL 
MAXVAL=VECTOR(  1 ) 
DO  10  J=2,MAX 

IF  (VECTOR(J)  .GT.  MAXVAL)  THEN 
MAXVAL  =  VECTOR(J) 
ENDIF 
10         CONTINUE 

END 

*    

SUBROUTINE  SMALL(TOTAL,VECTOR,LOW) 
* 

*  THIS  ROUTINE  TAKES  A  MATRIX  IN  VECTOR  FORM  AND 

*  FINDS  THE  SMALLEST 

*  VALUE  IN  EACH  ROW  AND  PUTS  IT  INTO  A  VECTOR 

INTEGER  TOTAL,ROWS,COLS 
REAL  VECTOR(TOTAL),LOW 
COMMON  /SML/  ROWS,COLS 

*  LOOK  AT  ONLY  ROWS 

LOW=VECTOR(l) 

* 

DO  10  1=1,  ROWS 

IF  (VECTORd)  XT.  LOW)  LOW=VECTOR(I) 

* 

10       CONTINUE 

END 
* 

SUBROUTINE  SHELL(MAXNUM,ARRAY,NO) 
* 

*  MAXNUM  IS  THE  NUMBER  IN  THE  VECTOR  TO  BE 


* 


SORTED 

INTEGER  MAXNUM,OFFSET,SW,NO,LIMIT 
DIMENSION  NO(MAXNUM) 
REAL  ARRAY(MAXNUM),TEMP,TEMP2 
LOGICAL  SWITCH 


* 


PRINT  *,  'SUB  SHELL'    OFFSET=MAXNUM/2 
1  IF  (OFFSET  .GT.  0)  THEN 

LIMIT=MAXNUM-OFFSET 

5  SWITCH=  .FALSE. 

* 

DO  15  J=1,LIMIT 
IF  (ARRAY(J)  .GE.  ARRAY(J+OFFSET))  THEN 
TEMP  =  ARRAY(J) 
TEMP2  =  NO(J) 
ARRAY(J)=ARRAY(J+OFFSET) 
NO(J)=NO(J+OFFSET) 
ARRAY(J+OFFSET)=TEMP 
NO(J+OFFSET)=TEMP2 
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SWITCH  =  .TRUE. 
SW=J 
ENDIF 
15  CONTINUE 

LIMIT  =  SW  -OFFSET 
IF(SWITCH)  GOTO  5 
OFFSET=OFFSET/2 
GOTO  1 

ENDIF 

* 

END 


SUBROUTINE  RAND(ASEED,TOTAL,RNDM) 

*  

*  Random  number  generator  -  Uniformly  distributed  in  (0,1) 

*  SEED  in  (1,247483647) 
DOUBLE  PRECISION  RANDM,ASEED 
INTEGER  TOTAL 

REAL  RNDM(TOTAL) 

* 

DO  20  1=1, TOTAL 
ASEED=ASEED+1.0 

RANDM  =  DMOD(16807.0D0*ASEED,2147483647.0D0) 
ASEED=RANDM 

RNDM(I)  =  SNGL(RANDM*4.6566128752458D-10)-1.0E-07 
20       CONTINUE 

ASEED=ANINT(ASEED) 

END 

*  _ 


SUBROUTINE  RANDN(ASEED,RNDM) 


*  Random  number  generator  -  Uniformly  distributed  in  (0,1) 

*  SEED  in  (1,247483647) 
DOUBLE  PRECISION  RANDM,ASEED 
INTEGER  TOTAL 

REAL  RNDM 

ASEED=ASEED+1.0 

RANDM  =  DMOD(16807.0D0*ASEED,2147483647.0D0) 

ASEED=RANDM 

RNDM  =  SNGURANDM*4.6566128752458D-10)-1.0E-07 

ASEED=ANINT(ASEED) 

END 
* 

SUBROUTINE  NORMAL(ALPHA,BETA,N,RAND2, 

+SLACK,MAXSLK) 

*    

*  PERFORMS  A  SIMULATION  OF  EXPERIMENTS  USING  THE  * 

*  WEIBULL  MODEL 

INTEGER  N 

REAL  ALPHA,BETA,MEAN,MAXSLK,RAND2(N), 
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+SLACK(N),RNDM 
MAXSLK  =  0.0 
PRINT  *,  'SUB  NORMAL' 
*--  CALCULATE  THE  STRENGTH  OF  THE  WEIBULL  FUNCTION, 
DO  10  1=1  ,N 
RNDM  =  RAND2CI) 

S  L  A  C  K  (  I  ) 
+EXP((LOG{-LOG<1-RNDM))+ALPHA*LOG<BETA)VALPHA) 
IF  (SLACKd)  .GT.  MAXSLK)  THEN 
MAXSLK  =  SLACK(I) 
ENDIF 
10         CONTINUE 

END 
*  END 


* 


SUBROUTINE 

+  UNIFM(RANGE1,RANGE2,N,RAND2,SLACK,MAXSLK) 
* 

INTEGER  N 

REAL  RANGE,MEAN,RAND2(N),SLACK(N),MAXSLK 

* 

CONST=RANGE2-RANGE  1 
MAXSLK=0.0 
DO  10  I=1,N 

SLACK(I)=RAND2(I)*CONST  +  RANGE1 
IF  (MAXSLK  .LT.  SLACK(D)  THEN 
MAXSLK=SLACK(I) 
ENDIF 
10       CONTINUE 
END 


SUBROUTINE  WEIBUL(ALPHA,BETA,TOTAL,BUNDLE,P) 


* 


*  PERFORMS  A  SIMULATION  OF  EXPERIMENTS  USING  THE 

*  WEIBULL  MODEL 

*  INITIALIZE 
REAL  ALPHA,BETA 
INTEGER  TOTAL,MAX 

REAL  P(TOTAL))BUNDLE(TOTAL),RNDM 
PRINT  *,  'SUB  WEIBUL' 
*--  CALCULATE  THE  STRENGTH  OF  THE  WEIBULL  FUNCTION, 
DO  10  1=1,  TOTAL 
RNDM  =  BUNDLE(I) 

P(I)  = 

+  EXP((LOGK-LOG<1-RNDM))+ALPHA*LOG{BETA))/ALPHA) 
10         CONTINUE 
END 


98 


LIST  OF  REFERENCES 


1.  Jones,    H.G.,    Mechanics    of    Composite    Materials,    pp.  1-28,     Hemisphere 
Publishing  Corp.,  1975. 

2.  Morely,  J.G.,  High-Performance  Fibre  Composites,  Academic  Press  Limited, 
1987. 

3.  Hoskin,    B.C.    and    Baker,    A. A.    (eds.),    Composite    Materials    for   Aircraft 
Structures,  American  Institute  of  Aeronautics  and  Astronautics,  Inc.,  1975. 

4.  Wu,    E.W.,    Class   Notes;   Reliability   in    Structures   and   Materials,    Naval 
Postgraduate  School,  Monterey  California,  1988. 

5.  Jones,   M.C.,   Composite  Reliability  via  Preloading,   Masters  Thesis,  Naval 
Postgraduate  School,  Monterey,  California,  September  1988. 

6.  Bell,    D.K.,    Composite    Reliability    Enhancement    via    Preloading,    Masters 
Thesis,  Naval  Postgraduate  School,  Monterey,  California,  June  1987. 


99 


BIBLIOGRAPHY 

1.  D'Agostino  R.B.  and  Stephen  MA.,  Goodness -of -Fit  Techniques,  Statistics: 
Textbooks  and  Monographs,  v.  68,  Marcel  Dekker,  Inc.,  1986. 

2.  Pachner,  J.,  Handbook   of  Numerical  Analysis  Applications,   McGraw-Hill 
Book  Comp.,  1984. 

3.  Thorndike,  R.M.,  Ph.D.,  Data  Collection  and  Analysis,  Gardner  Press,  Inc., 
1982. 

4.  Young,    H.D.,    Statistical    Treatment    of  Experimental   Data,    McGraw-Hill 
Book  Co.,  1962. 

5.  Composite   Materials:   Testing  and  Design    (Second   Conference),   American 
Society  for  Testing  and  Materials,  1971. 


100 


INITIAL  DISTRIBUTION  LIST 


1.  Commander 

Naval  Air  Systems  Command 

Assistant  Commander 

for  Sys  and  Engineering  (NAIR-05) 

2.  Defense  Technical  Information  Center 
Cameron  Station 

Alexandria,  Virginia  22304-6145 

3.  Library,  Code  0142 
Naval  Postgraduate  School 
Monterey,  California  93943-5000 

4.  Department  Chairman,  Code  67 
Department  of  Aeronautical  Engineering 
Naval  Postgraduate  School 
Monterey,  California  93943-5000 

5.  Dr.  E.M.  Wu,  Code  67Wt 
Department  of  Aeronautical  Engineering 
Naval  Postgraduate  School 
Monterey,  California  93943-5000 

6.  Lt.  Joseph  A.  Schmidt,  USNR 
4011  Chaney  Cove  Ct 
Dunkirk,  MD  20754 

7.  Dr.  Robert  Badaliance 

Naval  Research  Lab,  Code  6380 
Branch  Head/Mechanics  of  Materials 
Washington,  D.C.  20375 

8.  Dr.  Joseph  M.  Augl 

Naval  Surface  Weapons  Center 

White  Oak 

Silver  Spring,  Maryland  20910 


No. 

Copies 


101 


Thesis 

S337234  Schmidt 
cl 


Composite  reliability, 


Thesis 

S337234   Schmidt 

c.l        Composite  reliability, 


